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ABSTRACT 



> 

lO ' We use a high-temperature chemical network to derive the molecular abundances 

■ in axisymmetric accretion disk models around active galactic nuclei (AGNs) within 
100 pc using simple radial and vertical density and temperature distributions motivated 

^ by more detailed physical models. We explore the effects of X-ray irradiation and cos- 

■ mic ray ionization on the spatial distribution of the molecular abundances of CO, CN, 
CS, HCN, HCO"^, HC3N, C2H, and C-C3H2 using a variety of plausible disk structures. 
These simple models have molecular regions with a layer of X-ray dominated regions, 

^ ■ a midplane without the strong influence of X-rays, and a high- temperature region in 

c3 . the inner portion with moderate X-ray flux where families of polyynes (CnH2) and 

cyanopolyynes can be enhanced. For the high midplane density disks we explore, we 
find that cosmic rays produced by supemovae do not significantly affect the regions 
unless the star formation efficiency significantly exceeds that of the Milky Way. We 
highlight molecular abundance observations and ratios that may distinguish among 
theoretical models of the density distribution in AGN disks. Finally, we assess the 
importance of the shock crossing time and the accretion time relative to the forma- 
tion time for various chemical species. Vertical column densities are tabulated for a 
number of molecular species at both the characteristic shock crossing time and steady 
state. Although we do not attempt to fit any particular system or set of observations, 
we discuss our models and results in the context of the nearby AGN NGC 1068. 
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Introduction 



The growth of active galactic nuclei (AGNs) is an important key to the evolution of galaxies, 
and its connection with star formation inside galaxies remains of great interest in astrophysics. 

There are many unknowns in our knowledge of the accretion process and the structures of 
disks, such as mass accretion rates and disk thicknesses on scales ranging from fractions of a 
parsec to tens of parsecs. Other unknowns include the conditions for star formation in these en- 



eral models have been proposed to explain these quantities (e.g.. 


Shlosman et al. 1990; Goodman 


2003: Sirko & Goodman 


2003; Thompson et a 


[.boos : Chans etal. 


2007 


KroliklboOVilLevin 


2003; 


Wada & Norman 2007t 


VoUmeretal. 2008: 


Kawakatu & Wada 


2008 


. ISchartmann et al. 


20101; 


Hopkins & OuataerdboiO). 



The detection of molecules in extemal galaxies can help to characterize their star forma- 
tion regions and the physical conditions of the gas. Star formation takes place in a molecu- 
lar region as opposed to one that is mainly atomic or ionized. Carbon monoxide has been a 
main tool to obtain the molecular mass, whether in galactic clouds or extragalactic sources. In 
addition, many other molecular species have been detected in extragalactic sources, and show 
different features in different types of galaxies such as AGN-dominated g alaxies, (ultra- ) lumi- 
nous infrared galaxies ((U)LIRGs), and starburst galaxies. For example, iKohnd (120051) intro- 
duced the HCN/HC0'^(1 -0) line rati o as a diagnostic be t\ yeen AGN- and st arburst-dominated 
galaxies. Based on the observation by lTacconi et al.l (119941) . lUsero et al.l (12004) suggested a high 
HCN/CO intensity ratio in the AGN-containing galaxy NGC 1068. Unlike AGNs and starburst 
galaxies, (U)LIRGs such as Arp 2 20 and NGC 4418 seem to have more abundant complex organic 
species such as C9H9 a nd HC3N duahuis et allboOTl : kalto et allboOTl : Icostagliola & Aaltoll201ol : 
Costagliola et alboul) . 



The observed intensities are affected both by the abundances of molecular species and by 
the physical conditio ns, which in tu r n affe ct radiative transfer. The 3D-radiative transfer simula- 
tion of NGC 1068 by lYamada et aD J2OO7I) suggested that the higher line ratio of HCN/HCO+ in 
AGN-dominant galaxies is caused by the difference in chemical abundances rather than by exci- 
tation conditions. Chemical abundances are themselv es affected by different ph ysical parameters 
such as the ionization rate, temperature, and density. Lepp & Dalgamo ( 1996 ) examined the ef- 
fect of X-rays to explain the enhanced HCN/CO abundance ratio in NGC 1068 and the Galactic 
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Center. iMeijerink & SpaansI (|2005|) developed steady-state gas-phase chemical models of photon- 
dominated regions (PDRs) and X-ray dominated regions (XDRs), and w orked out an extensiv e 
grid of such model s to characterize extragalactic and galactic sources (|Meijerink et al.l |2007|) . 
Loenen et alj ( 20081) examined the effect of mechanical heating to enhance the HCN abundance 
in (U)LIRGs. IPapadopoulosI (120 lOh proposed an elevated cosmic -ray ionization rate as a possible 
dominant mec hanism in ULIRGs to cr eate cosmic -ray dominated regio ns (CRDRs), and its effect 
is discussed in IMeijerink et alj (|201l[) and IPapadopoulos et al.1 (1201 lah. T he spatial dependence 



of the effect of X-rays in NGC 1068 was studied by lOalliano et all (|2003|) . who modeled the H2 
ro-vibrational emission and the CO rotational emission using their X-ray irradiated model, which 
they find can explain the overall shape of the observed double-peaked emission. 

In this paper, we present detailed model calculations of the molecular abundances of a large 
number of species in the accretion disk of an AGN, using a black hole mass of Mbh = 10^ Mq, and 
assuming the gas density estimated from a Toomre-stable disk (see Section lZTI for details) out to ~ 
100 pc from the AGN core. This type of analysis should help us to use molecules as better probes 
of the physical conditions in AGN disks. In principle, when combined with chemical modeling and 
radiative transfer, molecular line observations can powerfully constrain the disk structure, density 
profile, and temperature profile on parsec to tens of parsecs scales in AGN disks, particularly 
when very high resolution images from ALMA become available. Unlike earlier authors, we 
focus our analysis on the radial and height dependence of molecular abu ndances, which to the best 
of our knowledge has never been done fo r AGN disks except for CO (IGalliano et al.ll2003[) and 
H2 (IGalliano et al.ll2003l : IWada et al.ll2009l : IPerez-Beaupuits et al.ll201l|) . The radial dependences 



of mo lecular abundances in the disk midplane were examined in our previous paper ([Harada et al 



20101) . Adopting parameterized vertical density structures, we add here the height dependence, an 
explicit calculation of X-ray and cosmic-ray ionization rates, and a self-consistent determination 
of the gas temperature including molecular line cooling. As a function of density and temperature 
throughout the disk we provide chemical abundances and average column densities for input to 
future radiative transfer calculations, and as a guide to interpreting current and future molecular 
line observations of AGN disks. 

The remainder of the paper is organized as follows. In Section [2l we discuss the range of 
physical conditions and parameters such as the temperature, density, and ionization rate used in 
our models, while in Section |3l we present our results. Section |4] contains a comparison between 
our results and observations, especially for molecules such as HCN, CN, and HCO"^, which have 
been observed at higher spatial resolution than other molecular species in NGC 1068. In addition, 
we include a discussion on how our results change as a function of different AGN disk models, 
and in different types of galaxies. We summarize our results in Section |5l 
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2. Model 

The chemical network used here is the OSU high-temperature network, which can be em- 
ployed up to a temperature of 800 K. With certain exceptions, such as the formation of H2 on 
granular surfaces, this is a gas-phase network. This network contains X-ray and cosmic -ray ion- 
ization and photodissociation from Lyman-Wemer photons and Lyman-a photons. Doubly ionized 
species are not included since we are interested in regions that are mostly molecular, where these 



Harada et al. 


(I2OIO) and 


Harada et al. 


(2012 



at pittp ://www .phy sic s . ohio- state .edu/^eric , and it will be updated on KIDA (Kinetics Database 
for Astrochemistry)Q, kida.uva.2012. As initial conditions, we used the so-called "low metal" el- 
emental abundances in which some depletion on grains is assumed for metals. There is also the 
possibility of higher abundances of heavier elements, approaching those of solar abundances. The 
initial abundances both in the low-metal and solar abundance cases are shown in Tabled] Physical 
conditions are estimated via methods described in the following sections. Important symbols are 
summarized in Table |2l 



2.1. Density Distributions 



As discus sed in our earlier paper ([Harada et al.ll2010|) . the gas density, p, at the midplane is 
assumed to be ([Thompson et al.ll2005|) 



^2 



where the Keplerian rotation frequency Q is given by 



(1) 



GM, 



(2) 



Here Mbh is the mass of the central black hole, a is the velocity dispersion of the stellar bulge, 
and r is the distance from the AGN core. We use Mbh = 1.5x10'Mq and a = 150 km s"' 
(|Greenhill & Gwinnlll997l: iTremaine et al.l l2002h . For a gaseous disk, the Toomre stability pa- 
rameter Q (|Toomrelll964|) is expressed as 



Q 



(3) 



' http : //kida . ob s . u-borde aux 1 . fr/models 
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where v., is the sound speed, n is the epicyclic frequency, and E is the surface density. For 2 < 1 
the disk is gravitationally unstable, while for Q> 1, the disk is stable. M any authors have proposed 



that t he value of Q should self-regulate to a value close to unity (e.g., |Paczynskilll978 



Gammie 



200 Observationally, Q does not dif f er from uni t y by m ore than a factor of a few over t he entire 



galactic disk according to lOuirkI (| 1 9721) . iKennicutti (| 1 9 89h . and lMartin & Kennicutt 



VoUmer et al. 



(2001). For the 



(120081) suggest 



accretion disk of NGC 1068, the vali ie of the Toomre pa ram eter that 

is ~ 5 based on the observations of iHicks et al.l (|2009l) and lOavies et al.l (I2007|) on 30 pc scales. 



There are also opposing views in the literature suggesting that Q can b e much higher than unity 
(e.g., Q ~ 100 in the Galactic Center: IvoUmer et al.lboosi) . Additionallv. lwada & NormanI hoOl\) 
questioned the use of the Q-parameter as a criterion at all. Nevertheless, we use two values of Q 
for our models, Q = I and 5. The midplane densities for 2=laren = 2xl0^ cm"^, 1x10^ cm"^, 
and 1 X 10"^ cm"^ at r= 1, 10, and 100 pc, respectively, while for 2 = 5, they are simply a factor of 
5 lower (see equation [T). 



p(r,z) = p(r,0)exp 



(4) 



The densi t y at h eight z from the midplane for a thin disk in hydrostatic equilibrium is given 

by lFrank et al. (Eooil to be 

2 \hJ ^ 

where h is the scale height, which can be defined by the relation h = S(r)/p(r,z = 0)j^ with S 
the surface density as a function of radius. The ratio of scale height to radius can be obtained 
([Thompson et al.ll2005|) from the expression 



r 23/2 ' 



(5) 



where /g is the fraction of the mass of gas over the total mass (gas + stellar). Since the gas fraction 
is uncertain, so is the disk scale height, and we varied h/r in our calculations, in developing four 
different disk models. These models, shown in Tabled are defined in three radial regions: < 1 pc, 
1 - 3 pc, and 3 - 100 pc. 



Model 2 is a relatively thin disk as in [Thompson et al.l (|2005h . in which equation dH) is used 
for the dependence of density on height. In this model, the value ofh/r first increases quadratically 
with radius then increases linearly. This model is motivated in part by the recent h y drody namic 
simulations of star-forming disks with a central black hole by iHo pkins & Ouataert ( 2010 ). who 
show h/r m their models as a function of various parameters, including the gas fraction, surface 



- The vertical density distribution adopted in equatio n dUi is only a rough approximation, and self-consistent m odels 
may yield a different functional form for p{z) (see e.g., Thompson et al. 2005 : KrohklEoO? : Chang et al ] l200i . 

^For Model 4 and the torus portion in Model 1, discussed later in this section, h is defined in a different way, and 
the surface density is lower than in other models. 
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density, and the central black hole mass. Figure 10 in their paper shows that /z/r < 0. 1 on 1 - 100 pc 
scales in most cases, although there are fairly large fluctuations. 



The other three models we consider contain tori of assorted radial widths. IJaffe et al.l (|2004|) 
observed a torus of /z/r ~ 0.6 in NGC 1068 on a few pc scales, and our Model 1 is a thin disk model 
analogous to Model 2 but with a torus of h /r = 0.5 located at 1 - 3 pc. A relatively thick torus has 
also been suggested by iHicks et al.l (|2009h . with h/r > 1 on 30 pc scales based on their studies of 
ro-vibrational H2 emission. This warm H2 traces molecular gas at a few thousand Kelvin, and it 
is not certain if the warm gas is tightly coupled to the colder midplane of the disk. Nonetheless, 
we also consider geometrically thick torus cases in Models 3 and 4, with h/r = 0.5 for radii of 
1 - 100 pc. We use the height dependence expressed in eq. (H) for Model 3. However, this is a 
thin disk approximation, and may not express the correct height dependence of density in a thick 
torus, which could also be affected by t he irradiation from t he AGN core. We therefore employed 
the height dependence of density from IChang et al.l (|2007|) by fitting their numerical expression 



to p(z) = p(0) for z < 0.02/z, p(z) oc r"^-^^ for 0.02/z < z < h, and piz) oc exp(-^) for z> h (see 
also footnote 2). Fig [T] shows the densities of all the models plotted with contours representing 
column densities from the X-ray source. The total molecular mass is ~ 7 x 10'Mq for Model 1 
and 2, ~ 3 X 10^ Mq for Model 4 within 30 pc. The total m asses of Model 1 and 2 are larger 
than those in the hydrodynamic models of IWada et al.l (|2009|) . where the total molecular mass is 
~ 1.5 X lO^Mgwhen the black hole mass is 1.3xlO''M0. 



2.2. Temperature 

We calculate the gas temperature from the heating and cooling rate; the coupling with the 
dust temperature is also included. For the range of dust temperatures in the AGN disk, we used the 
blackbody approximation, determined by the luminosity of the central AGN: 

-1/2 / J ^ ^1/4 

(6) 



\pc / \2 X 10^^ erg s ^ J 



where the normalizing luminosity of 2x 10"^^ erg s ' derives from the observation by iMason et al. 



(120061) of NGC 1068. Using this formula, T^ust = 750, 240, and 75 K at r = 1, 10, and 100 pc, 
respectively. To calculate the gas temperature, we include heating from both X-rays and cosmic 
rays. Dust-gas collisions provide gas heating or cooling, depending on the dust-gas temperature 
ratio. The cooling rate is determined by molecular line cooling. Although photons produced by 
molecular lines can in principle escape in any direction, we make the simplifying assumption that 
the cooling line photons always escape from the AGN disk vertically, which reduces the prob- 
lem to a plane-parallel geometry at every radius. Even so, the calculation of the gas temperature 
is complicated by the fact that the chemical abundances, which determine the cooling rate, and 
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temperature are interdependent. Our procedure is as follows. We first calculate the chemical abun- 
dances assuming that the gas temperature is T{r) = T^ustif) everywhere. Given these abundances, 
we then calculate the temperature in thermal equilibrium, with the heating and cooling rates de- 
scribed in Appendix |A] Once the temperature is determined, from this initial step, the abundances 
and cooling rates are updated, and then the temperature is recomputed. This procedure is iterated 
until convergence. 

Figure |2] shows the gas temperature at each part of the disk in Models 1-4. In the regions 
that are irradiated strongly by X-rays, the temperature can reach T > lO^K. For columns of 
~ 10^"* cm"^. X-rays are attenuated, and the gas is predominantly heated by other sources. When 
r < 5pc, the gas temperature is coupled with the dust temperature because of the high density. 
Farther away from the AGN core, the temperature slowly decouples from the dust temperature and 
goes down to lOK at r = lOOpc. We calculated molecular abundances and followed the iterative 
procedure described above only wh en T < 5000 K since few molecules exist at higher tempera- 
tures. Observations of ammonia by lAo et al.l (|201l|) suggest warm ~ 80K and ~ 140 K molecular 
gas components on a 1.2 kpc scale, which implies that other sources of heating such as turbu- 
lence/shock heating may be dominant on these scales, as has been suggested by many authors. 



Our temperature calculation was checked against the results by IMeijerink & SpaansI (|2005|) . 
We repr oduce the general behavior of the run of temperature in XDRs as a function of column, 
although iMeijerink & SpaansI (|2005|) use a more detailed calculation of the X-ray penetration cross 
section than employed here, which makes X-rays in our results attenuate somewhat faster. 



2.3. Timescales 



Although the density and temperature of our disk models are independent of time, several 
timescales, and their ratios, determine in part the abundances that might be observed in AGN disks. 
The relevant timescales are the accretion time ?acc, the shock crossing time Across, and the chemistry 
time ^chem- The accretion time, tacc, is the time for the medium at radius r to be accreted. In viscous 
disks, ?acc ~ (al7)~Hr//i)^ where a is the visco s ity pa rameter. With a ~ 0.01 -0.1, ?acc ^1000 
for /i/r ~ 0.1. However, iHopkins & OuataertI (|2010l) show that the gravitational torque from the 
stars on the gas is likely to be the dominant angular momentum transport mechanism in disks that 
feed AGN. In such disks tacc ~ where \a\ is the fractional magnitude of the asymmetry 

dLvnden-Bell & Kalnaisl[l972l) . In the simulation bv lHopkins & OuataertI hoi(h . \a\ ~ 0.01 -0.3, 



implying that ?acc ~ 3 - 10017 ^ 

The accretion time in any model should be compared \yith the shock crossing time, fcross- 
The medium is turbulent, and the simulation by IStone et al.l (|l998|) shows that the dissipation 
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time of shock waves should be on the order of VT^ . The velocity of these shocks can be tens 
of km s~\ and they can disso ciate the molecules, effectively resetting the molecular abundances 
(|Hollenbach & McKedllQSQ ). For the disk models we consider, the shock crossing time Across = 
5 X 10"^ yr xQ~'/^( ^q^"^_3 )~'^^ and the values in the model are 5x 10^, 4x 10"^, and 3x 10^ yr at r = 
100, 10, and 1 pc. 

There are three cases possible for the relationship of these three time scales: (1) ?chem < Across < 
tacc, (2) Across < ^chem < ?acc, and (3) Across < kcc < ^chem- If the shock wavcs are dissociativc, then 
starting from our initial abundances at ? = 0, the chemical abundances calculated at t = Across will 
represent average time-independent abundances. Conversely, if the shocks are not dissociative and 
molecules survive shock passage, then the chemistry calculated to t = ?chem (i-C-, the true steady- 
state time) should be representative of the abundances as long as t^cc > ^chem- This assumption is 
valid for cases 1 and 2, but the change in physical conditions due to accretion should be considered 
for case sj^l 

Starting from our initial abundances at f = 0, choosing the so-called low-metal values, (see 
Table[il), we calculate the time evolution of the abundances at each point in density and temperature 
in our disk models until both Across and ?chem- Depending on the value of ?acc/^chem for the molecule 
in question, the effect of accretion is sometimes important for smal l enough r. We disc uss these 
cases explicitly. The accretion in the midplane was considered in iHarada et al. ( 2010 ). and its 
effect is discussed in later sectio n s. As discussed in models of protopl anetary disks (see e.g., 
Morfill& Voel5ll984l : lOaillboOll: lllgner et all [2OO4I: Iwillacv et all 12006 ). the effects of vertical 
and radial mixing in addition to accretion might be important, but consideration of these factors is 
left to future work. 



2.4. Ionization Rate from Cosmic Rays 

Primary cosmic ray protons are an important source of ionization in star-forming disks, and 
their secondary particles such as photons, gamma rays, electrons and positrons can also cause 
ionizations (primary electrons are less energetically important than primary protons and their 
secondaries, and thus these a re neglected here, as in standard analogous Galactic calculations; 



Umebayashi & Nakanolll981l) . Cosmic rays are likely produced by the supernova explosions and 
stellar winds of massive OB stars, and thus their production rate is directly connected with the star 
formation rate. Cosmic rays can also come from the central AGN. 



'^We are aware that these are two extreme cases. In between there is a C-shock, which dissociates only some of the 
molecules. 
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2.4.1. Volumetric Star Formation Rate 



The star formation rate per unit volume, p'*, may be written as (|Kruniholz & McKedl2005l) 

P _ ..^3/2(-l/2 



^dyn 

2 X 10"^ 



Mo 



pc3 yr VO.Ol/ VlO^cm 



3/2 



(7) 



where ?dyn ~ 1/ \/Gp and is the star format ion efficiency, whose value of ~ 0.01 in the above is 
suggested by the results of iKennicuttI (Il998[) . Since the star formation efficiency and cosmic ray 
production rate in A GN disks is iincertai n, we vary y from 10~^ - 10~^. The star formation rate per 
unit area inferred by loavies et all (|2007|) within 35 pc from the AGN core in NGC 1068 is IOOMq 
yr"^ kpc"^. For the disk density models we use (see eq. [T]), a star formation efficiency of 10"^ - 10"^ 
corresponds to this value of the star formation ra te per unit area. Such low values for u at small 
radii may be consistent with models of feedback ([Thompson et al.ll2005l : Davies et al.ll2007n . 



2.4.2. Cosmic Rays from Supernovae & Stellar Winds 



In the Milky Way, the supern ova rate is a few per centurv (iDiehl et al. and the estimated 

star fo rmation rate is ~ 1 yr"' (lGilmorell200ll: iMcKee & Williamslll997l : lRobitaille & Whitney 
2010h . Assuming that every 100 Mq of star formation produces 1 supernova, and that a fraction 
/cR of the supernova explosion energy (Esn = 10^' erg) is injected into the ISM per explosion, the 
cosmic ray luminosity per volume (Lcr/V) is directly proportional to (eq.|7]). The energy density 



of co smic rays, which is directly connected with the cosmic ray ionization rate (|Papadopoulos 
2010l) . can be written roughly as 

LcR 



V 



■ min[?escape : ^pp] 



(8) 



where ^escape is the escape timescale due to either diffusive (as in the Galaxy) or advective losses 
(as perhaps in starburst winds), and tpp is the timescale for inelastic proton-proton collisions (pion 
production). If escape is dominated by an outflow, then 



«' escape 



■ ^wind 



h 

V 



3 X 10^ yr 



(^) 



/300km/s 



(9) 



wher e h is the gas scale height and v is an assumed wind velocity. This should be compared with 
(e.g.. ISchlickeisei]l2002[ lTorresll2004l) 



tpp ~ 



/lO^cm-3 
50 yr 

V « 



(10) 
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Although only a rough estimate, the fact that tpp <^ ^wind suggests that in the disks we consider 
cosmic ray losses are strongly dominated by pion production losses. Hence, 



UcR - fcR P* fsN EsN tpp 

~ 4 X 10"^ erg cm"^ 



lO^cm- 



1/2 



0.01 



(11) 

(12) 



where we have taken /sn = 1 SN/IOOMq, E'sn = 10^^ e rg, and fcK = 0.1. Taki ng the cosmic -ray 
energy density of our Galaxy to be ~ 2 x 10"^^ erg cm"^ ([Boulares & Coxlll990() . 



U, 



CR 



2x 10^ 



CRfial 



lO^cm- 



1/2 



0.01 



(13) 



Althoug h simplified^ this res ult for Uqk as a function of n is in accord with the more detailed mod 
elling of Lacki et al. ( 2010al) in high-d ensity starburst galaxies (their Fig. 15). See Papadopoulos 
(|2010() and lPapadopoulos et al.l (|2011br) for a more detailed discussion of the effects of high cosmic 
ray ionization rates on the interstellar chemistry and thermodynamics in dense starburst galaxies. 

Note that with the average density of the Milky Way (hmw) ~ 1 cm"^, equation (fT3]) implies 
2. The value from this simple approximation is expected to be higher than the Galactic 



UcR.Gal 



value since the majority of the cosmic rays escape from the Galaxy via diffusion before energetic 
losses (e.g.. IStrong et al.ll2000|) . Taking the galactic value of the cosmic-ray ionization rate to be 
1 X 10"^' s"\ we have that 



14^-1 



CcR = 2xlO-^^ 



106 



cm 



1/2 



0.01 



(14) 



where we have again assumed /cr = 0. 1 and Esn = 10 erg- The value of the cosmic-ray ionization 



rate in the midplane with u = 0.01 is Ccr 
linearly scales as u. 



10"^^ s"^ at r ■■ 



1 pc , and 10 



"'^ s ^ at r : 



lOOpc, and it 



Stellar winds from massive OB stars can also be a source of cosmic rays in the dense ISM of 
AGN disks. The power input in the Milky Way is thought to be a bout an order of magnitude lower 
than t hat supplied by supemovae, and so we neglect them here dCasse & Paulll 19821 : ISchlickeisei 
20021) . However, even if they were to be a significant source of cosmic rays comparable to su- 
pemovae, this would only contribute at order unity to the ionization rate, and result in an overall 
rescaling of the constant in equation (fT4l) since the supernova rate and massive OB star birth rate 
scale directly with the star formation rate. 
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2.4.3. Cosmic Rays in the AGN core 



The A GN core is another p ossible source of cosmic rays. Assuming that the gamma-ray flux 
reported bv Lenain et al. ( 2010 ) is hadronic (from vr" decay), the ratio of the gamma-ray flux to the 
bolometric flux i ndicates that NGC 1068 possesses near the maximum attainable ratio from star 
formation alone (ILacki et al.ll2010bh . Thus, a dominant contribution to the total cosmic ray budget 
from the AGN core of 1068 is neither supported or excluded. 

Finally, even without the contribution to cosmic ray production from supemovae or OB star- 
formation, there should be some minimum flux of cosmic-rays produced b y radioactive nuclei. For 
a mini mum value, we use CcR,min = 10"^** s"^ based on ^''Al, as discussed in lUmebayashi & Nakano 
(I2OO9I) . 



2.5. Effects of stellar UV-photons 

When gas is irradiated by a nearby OB star, the UV-photons can dissociate molecules. If the 
star formation rate is high, the UV-photons can affect the chemistry significantly. At the same time, 
the penetration depth of UV-photons is much less than that of X-rays or cosmic rays; thus, its effect 
can be minimal in a dense medium. We use a simple approximation to determine if UV photons are 
important or not. The rate coefficient of photodissociation by UV-photons is k^i oc Goexp(-7Av) 
where 7 ~ 2-3, which means that a region with high visual extinction surrounding the star is only 
affected significantly quite near to the star, where we refer to the conditions as those of a photon- 
dominated region, or PDR. The ratio of the volume of a PDR where Ay < 10 to the entire volume 
(the inverse of the OB stellar density) can be estimated by 

VpDR _ . 47r(A^4,io)% 
^ - P*^;^foBroB 

where NavIO is the PDR column density, for which we use a value of 2 x 10^^ cm"^. When n = 
10^ cm"^, the PDR volume only accounts for a fraction of 2x10"^ of the entire region. If the 
star formation efficiency is higher, PDRs may affect the observable molecular abundances more 
significantly, but here we choose not to include UV photons in our model. We also ignore the UV 
photons from the AGN core because of the high obscuration. 
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2.6. X-rays from the AGN Core 



The ionization rate caused by X-rays interacting with H2 is given by lMaloney et alj (|l996l) via 
the equation 



,ec , a{E)F{E)dE 



(16) 



where the flux per unit energy, F{E) = FQ{E)e~'^^^\ N^ec is the number of secondary ionization 
events caused by electrons produced by the primary ionization, t(E) is the optical depth, and a(E) 
(cm^) is the primary ionization cross section at energy E, which is given by 



a(E) = 2.6 X IQ-^^E-^/^cm^ (1 keV < E < 7 keV), 
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a(E) = 4Ax IQ-^^E 



/^cm^ (7keV<E). 



(17) 



(18) 



We assume that Fo(E) oc E'^-'' (i Beckmann et al.1 120061) . and we normalize the flux at each energy 
from the total X-ray flux. The penetration of X-rays can be divided into two ranges of energy, 
above and below 1 1 keV. For energies lower than 1 1 keV, photoionization, whose cross section is 
given by equations (I17H18I ), is the dominant process limiting penetration. Note that although our 
gas-phase metal abundances are very low because of condensation onto grains, we still assume 
the total metal abundance is in the gas phase for the purposes of this calculati on, as if the metals 
are not locked on grain surfaces. As noted in IMorrison & McCammonI (119831) . the cross sections 
do not differ significantly even if m ost of the heavi er elements are condensed in grains, and thus 
the estimate of the cross section by iMaloney et al. I il996) above should still be valid. However, 
for higher metallicity than we employ here this approximation fails and the photoionization cross 
section should be computed more fully. 

For energies higher than 1 1 keV, Thomson/Compton scattering cross section dominates. Al- 
though the ionization rate is still dominated by photoionization for a column A^h ^ 10^"* cm"^ 
(|Maloney et al.ll 19961) . our models have higher midplane column densities. Since the differential 
cross section of Compton scattering has angular dependence, a precise calculation of the penetra- 
tion would require a more complicated transport program. In our calculation, we simply assume 
that the optical depth for high-energy X-rays is reduced by a factor determined via the assumptions 
of isotropic scattering and random walk. Until an energy of ~ 100 keV, the cross section given by 
Klein-Nishina formula does not differ much from Thomson cross section of cr 7- = 6.7 x 10"^^ cm^. 
Once the high-energy X-rays are reduced to 1 1 keV energy, the normal one-dimensional treatment 
is used. 

The total X-ray luminosity from observations is not the best measure of the true luminosity 
because of obscuratio n. The intrinsic X-ray luminosity for the A GN core of NGC 1068 is estimated 
to be lO"^^""^"^ erg s"^ (llwasawa et al.lll997l : IColbert et al.ll2002|) . In our calculations, we use total 
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luminosity of 6 x 10'^'' erg s ^ over 1 keV to 100 keV with a photon number flux with the energy 
dependence of E"^ °. 



2.7. The Relative Ionizing Effects of X-rays and Cosmic Rays 

For most regions of the disk, the ionization caused by X-rays exceeds that due to cosmic rays 
unless v is as large as 0.01 (see eq. [I4l). We therefore first discuss the disk as an X-ray dominated 
region (XDR). The ionization due to X-rays depends upon position and the density distribution of 
the disk model chosen. Figure |3] shows the ionization rate throughout our disk models caused by 
X-rays. Rather than plotting the ionization rate explicitly, we have chosen to plot C,x/n, since this 
parameter determines the steady-state chemical ab undances prov ided that direct dissociat ion by 



UV-photons is negligible (ILepp & Dalgamolll996|) . For example. Lepp & Dalgamd (|l996|) show 



that the fractional abundances of CN and HCN can be enhanced with moderate values of (x/n. In, 
particular, their peak fractional abundances are achieved when (x/n ^ several x 10"'^ cm^ s"^ At 
the lower value of Cx/>^ = 10"^' cm^ s~\ their abundances decrease by two orders of magnitude. 
When Cx/n > 10"'*^cm^s"\ the fractional abundances decrease with increasing C,x/n since the 
molecules are dissociated rapidly. 

The upper left and right panels in Figure |3]show C,x/n with 2=1 for Models 1 and 2, re- 
spectively. These are high-density thin-disk models with and without inner tori. Here it can be 
seen that the ionization caused by X-rays has trouble penetrating to the midplane. In these mod- 
els. X-rays penetrate significantly only into a layer on the surface of the disk and create a region 
where C,x/n ~ 10"^' - 10"'^ cm^ s"' that has typical XDR chemistry, but C,x/n < 10"^^ cm^ s"' near 
the midplane. In the lower left panel in Figure |3j on the other hand, higher ionization levels are 
found; the results are for a lower-density {Q = 5) version of Model 3, which contains a torus out 
through 100 pc. Although Cx/" < 10"^^ cm^ s"' around the midplane, there is a thick layer with 
Cx/n= 10-1^-10-21 cm3 s-\ 

Model 4, shown in the lower right panel, has a higher X-ray ionization than Model 3 since 
the density falls off faster with height above the disk. Since X-rays penetrate more because of the 
lower density, there are more regions with Cx/^ > lO-'^cm^ s~^, where molecules cannot exist. 
Since Models 1 and 2 have a higher h/r ratio at larger r, there is a more extended XDR at larger 
r for these models. If we compare Models 3 and 4, which have constant h/r for r > 1 pc, we find 
that Model 3 has a more extended XDR at smaller r, while X-rays can reach larger values of r in 
Model 4. 



For those regions where X-ray ionization is relatively small, such as the midplane of the disk 
at large radii, the ionization rate from cosmic rays can possibly dominate if v is relatively large 
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(see eq. [141) . For example, when 2=1 and v = 0.01, Cx/n ranges from 8 x 10"^^ cm^ s"^ at 1 pc to 
8 X 10"'^ cm"* s"^ at 100 pc from the AGN core. The effect of cosmic rays becomes more important 
for larger r because of the relatively lower density (see Section |3^ below). 



3. Results 

3.1. Models with X-ray Ionization Only 

In this section we explore the effects of X-rays without the inclusion of cosmic-ray ionization. 
The degree to which X-rays can affect the disk chemistry depends on their penetration into the disk, 
which depends on the disk structure. Far from the midplane, molecular abundances are very low 
given the very large X-ray ionization and lower density. Only regions with considerable molecular 
abundances are shown in the figures to be discussed below. In between the "no-molecule zone" 
and the midplane, there is an XDR, the extent of which is determined by the ability of X-rays to 
penetrate. At closer distances and higher temperatures, the temperature affects the chemistry sig- 
nificantly over a range of X-ray ionization rates, and we refer to these portions of the disk as "high 
temperature-synthesis" regions. The main results of this high-temperature-synthesis effect are to 
increase the abundances of some species such as HCN, the family of CnH2, and cyanopolyynes 
(e.g., HC3N). The effect is significant for T > 300 K at early times before CO fully forms, and can 
be seen in the results presented at the shock crossing time t = Across (see Section [23l) . Yet, in the ab- 
sence of X-ray heating, this effect still remains at steady- state after CO fully forms for T > 500 K. 
Importantly, though, in our calculations the regions where T > 500 K are in general strongly heated 
by X-rays, which tend to dissociate molecules, and therefore the high-temperature-synthesis effect 
on the final abundances is miti gated in regions of re latively low columns (< 10^"* cm"^). Details of 



this chemistry are discussed in lHarada et al.l (|2010|) . 



We find that the results for the abundances in Models 1 and 2 are not qualitatively different 
except for the values at 1 pc. However, Model 4 is quite different from the others. For these 
reasons, we discuss, compare, and contrast Model 1 with 2=1 and Model 4 with 2 = 5, with 
an emphasis on fractional abundances X of important molecules with respect to the total hydrogen 
abundance. As discussed in Section [2?T1 Model 1 is a high-density thin disk model in which X-rays 
have more difficulty penetrating than in Model 4, which is a relatively low density geometrically 
thick torus. 
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3.1.1. Results for Model 1 with Q=\ 

The fractional abundances of selected molecules are shown in Figure |4] at t = Across, and Figure 
[5]for steady-state. The fractional abundance of each molecule is discussed below. 

HCN : The fractional abundance of HCN att = Across and at steady state is shown in the upper 
left panel of Figures |4] and |5l respectively. As mentioned above, we started all elements except for 
H2 from atomic form as in Table [T] At t = Across, there is a fraction of carbon that is not yet locked 
up in CO. This free carbon helps to form a high abundance X(HCN) > 10"' when T > 100 K in 
most regions of the disk. 

Although the HCN abundance can be enhanced at early times in the calculation (before steady 
state is reached) when T > lOOK, the temperature in the midplane of AGN disk is in general not 
high enough (T < 500 K) to keep the HCN abundance high as steady state is approached except 
for r < few pc where the dust temperature is highest. When accretion is considered, there might be 
less of a high-temperature-synthesis effect observed in the HCN abundance, because if accretion 
is rapid (e.g., ~ 10- lOOQ"') ?acc can be less than fchem, and one does not expect the steady-state 
abundances to be reached. 

CN: The time at which the peak in CN abundance is reached is much earlier than that of 
HCN. Since HCN is more stable, CN eventually becomes HCN either through an indirect process 
involving ions at cold tem peratures or via a rn ore efficient direct reaction with H2 at high tempera- 
ture (see equation (21) in iHarada et al.lboiol) . X(CN) peaks at ~ 3 x 10"*^ on a timescale 10- 100 
times shorter than Across- By ? = fcrossj ^(CN) has decreased by orders of magnitude, and is only 
enhanced in the XDR layer where X(CN) ~ a few x 10"*^. In the high temperature regions, the CN 
abundance is lower because of the efficient production of HCN from CN. Within the XDR layer, 
the CN fractional abundance is much higher in the colder, more distant, part. On longer timescales, 
the steady-state value of CN near the midplane is orders of magnitude lower than the abundance at 
t = Across, but a similar degree of enhancement is still seen in the XDR layer (see Figures Hand |5l 
upper right panel). 

HCO^ : The fractional abundance of HCO"^, depicted in Figures |4] and |5l does not show 
strong time dependence. Once CO is formed, HCO"^ can be made by receiving a proton from a 
protonated ion such as H3. Its abundance increases in moderately ionized regions, and we find a 
higher abundance of HCO"^ in the XDR layer, similar in location to where the abundance of CN 
is enhanced. At both t = ?cross and steady state, the highest abundance is seen at inner disk radii 
far from the midplane and in the outer XDR region, where the fractional abundance of HCO"^ is 
larger than that of CN with a value of X(HCO+) ~ (3-10) x 10"**. Starting the simulation with 
solar elemental abundances, we find that X(HCO"^) is reduced, where X(HCO'^) < 10"^ even in the 
XDR layer, since atomic ions such as S"^, Na"^, and Mg"^ hold most of the positive charges instead 
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ofHCO+. 

HC3N, C2H2, C-C3H2, C2H: Similar to HCN, there is a high-temperature enhancement for 
HC3N and the family of C„H2, although the C2H abundance decreases with increasing temperature. 
As seen in Figure ID the fractional abundance of HC3N exceeds 10"^ when r < 5 pc at f = fcross- 
In steady-state, however, there is high-temperature-synthesis effect only at r < a few pc where the 
disk is quite hot, as can be seen in Figure |5l The fractional abundance of C2H2 at t = fcross also 
peaks at high temperatures at a value of about 10"^ for r ~ lOpc, and the steady-state abundances 
decrease by orders of magnitude as they do in HC3N. The abundance of C-C3H2 follows a similar 
trend as C2H2. 

The fractional abundance of C2H is orders of magnitude lower than that of C2H2 in the high- 
temperature-synthesis regions. For r ~ 100 pc, the fractional abundance of C2H is around 1 x 
10"^ at t = ?crossj and decreases by orders of magnitude on longer timescales as steady-state is 
reached. When calculated using solar elemental abundances, the abundances of these carbon-chain 
molecules peak at slightly higher values in high-temperature- synthesis regions simply because of 
the higher elemental carbon abundance, but are orders of magnitude lower for other regions. 

CS : In general, the CS abundance peaks at times before t = t^oss, after which most of the sulfur 
goes into SO or SO2 in an oxygen-rich (C/0< 1) environment. By the time t = Across ^ ^(CS) < 10"^ 
for almost all regions, and SO/CS>10 for r < 20pc. For a larger r, due to the lower density and 
lower temperature, the chemistry evolves more slowly, yet there is still a high abundance of CS. 

QH^ ; Since OH"^ can survive in regions that have higher X-ray ionization rates than other 
molecules, it is abundant in the XDR layers over a wide range of r. In Model 1, there is only a 
small region where 0H+ is enhanced because of the adopted density structure and the attenuation 
of X-rays from the torus. 



3.1.2. Results for Model 4 with Q = 5 

The biggest difference between Models 1 and 4 is the larger fraction of the XDR layer over 
the entire disk in the latter. Its presence changes the abundance pattern of a number of species. As 
depicted in Figure |6] an enhancement of CN compared with Model 1 (see Figure H]) corresponds 
generally to the expanded location of the XDR layer. Here, the CN abundance at t = fcross is ~ 
(1 -5) X 10"^, which slightly exceeds the steady-state value of ~ (1 -2) x 10"^. Other abundances 
and abundance ratios are affected by the lowered densities in Model 4. For r > 5 pc, the CS/SO 
abundance ratio becomes higher than unity at t = faoss- The enhancement does not correspond to 
the XDR region, and it is instead due to the slow chemical evolution because of the relatively low 
density. Finally, the additional X-ray flux can dissociate the carbon-chain molecules in regions that 
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correspond to high-temperature-synthesis regions in Model 1 . These species only have moderate 
abundances at f = tcmss, and their abundances decrease by orders of magnitude at steady state. 



3.1.3. Calculated Column Densities and Ratios 

Our models lead to the prediction of vertical column densities by suitable integration of the 
calculated abundances. We list the calculated vertical columns at six radii ranging from 1 - 100 pc 
for 22 atomic and molecular species as a function of radius and temperature. Column densities at 
t = ?cross are listed in Tables and |7]for Model 1 (2 = 1), Model 2(Q=\), Model 3(Q = 5), 
and Model 4 (Q = 5), respectively, while Tables [8l |9l [lOl and \TT\ contain results for the column 
densities at steady state. 

The calculated column densities are highly dependent on the disk models. The molecular 
column densities of Model 4 are much lower than in Models 1, 2, and 3 partly because of the larger 
"no-molecule zone" due to high penetration of X-rays, and also because the z-dependence of the 
density gives a lower total column density. For Model 1 at 1 pc. Model 3 at r < 3 pc, and Model 4 
at 1 pc, the molecular column densities are smaller than the column densities at larger radii because 
the molecules are dissociated by X-ray irradiation. 

From the calculated columns, we can obtain ratios relevant to observations for each of our 
models. Below we discuss the specific ratios CN/HCN, CN/CO, HCN/CO, and HCO+/CO. Al- 
though strictly ratios of column densities, we often refer to them as abundance ratios below. The 
ratio of vertical column density at t = Across and at steady-state is plotted vs r for CN/HCN in Fig- 
ure IT] In general, this ratio decreases with decreasing r due to the higher temperatures at smaller 
radii, which favor HCN over CN except for the very inner part where the regions contain very few 
molecules. These regions do not contribute to the total molecular column densities. Because of 
their larger XDR layers. Models 3 and 4 have the highest ratio of CN/HCN, approaching unity 
at 100 pc. The column density ratio for CN/CO is the highest in Model 4 because of the higher 
portion of XDR layers, ranging from lO"'*- lO""* at f = t cross and 10"^ - 10""^ at steady state for most 
of the disk. There is no strong radial dependence for r > 3pc (Fig. [8]). 

Figure |9] shows the HCN/CO column ratio att = Across and at steady state. The abundance ratio 
directly reflects the high-temperature enhancement of HCN att = ta-oss, except for the inner regions 
of Model 4, where the strong X-ray dissociation causes HCN/CO to decrease. The HCN/CO ratio 
is 10"^-0.01 for r < 10 pc, and decreases for larger r because of the decreasing temperature to 
several x 10""* for t = Across- At steady state, the HCN/CO ratio drops down to ^lO""^. As shown 
in Figure \T0\ the HCOVCO ratio is around 10"^- 10^ for most of the disk. Finally, Figure [TT] 
shows the ratio of HCO+/HCN. When the HCO"^/HCN ratio is too low, it is influenced by the high 
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HCN/CO abundance ratio. 



3.2. Effects of Inclusion of Cosmic Ray Ionization Rate 

Our results presented so far have been obtained without the inclusion of cosmic ray ionization. 
We have also run models that include a range of cosmic ray ionization rates by varying the star 
formation efficiency v (see Section The effects of cosmic rays from star formation become 
dominant over X-ray ionization only when v > 0.01. Figure [l2] shows the variation of the total 
vertical column densities of HCN, CN, and HCO"^ respectively as a function of radius for values 
of v equal to 0, 10""*, 10"\ and 0.01 at t = Across for Model 1, Q= 1. The HCO"^ abundance can 
be enhanced by a factor of a few when the star formation efficiency is higher, although HCN 
decreases when cosmic-ray ionization is included. The column of CN is enhanced only slightly 
when u = 0.01. 

Inclusion of cosmic-rays can also dissociate complex molecules such as C3H2 and HC3N for 
u > 10"^ (see eq.|7]). Although cosmic-rays can heat the gas to a few x 100 K, the high ionization 
rate can cause dissociation of molecules preventing more complex molecules from forming. 



3.3. Line intensities: Preliminary Results 

We have computed observable intensities o f selected molecules using th e publicly available 



three-dimentional radiative transfer code LIME (iBrinch & Hogerheijdell2010h . Since our models 



span a wide range of physical conditions, including those that are highly optically thick, we are 
unable to achieve complete convergence in some regions, for some species and transitions. In 
particular, in LIME the signal-to-noise ratio of the population convergence is a ratio of the current 
population over the geometric mean of the fluctuations in the population over the last 5 iterations. 
Out of 80,000 grid points used, the lowest the signal-to-noise ratio is below 1. The average value 
of signal-to-noise ratio for CO reaches more than 100 although these ratios are lower for HCN and 
HCO"^, which is somewhere above 40 for model 2 and 20 for model 4. Lower levels of convergence 
come mostly from less populated high J levels with weaker emission, and should not affect the 
overall results for the low-7 lines (< 4) discussed in this paper. Even so, the results presented here 
are thus preliminary, and should only be used to understand general trends. With this precaution, 
we present our preliminary results on the line transfer here with selected models att = tcmss- 

The LIME code requires density, temperature, molecular abundance, velocity, and the doppler 
width for the input parameters. The velocity of the rotation and doppler width of hil are used and 
images are produced with a viewing angle of 45 degrees. Results presented here are convolved with 
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a bea msize of 0.25 arcsecond (~ 18pc at the distance of NGC 1068) using MIRIAD (|Sault et al 



1995 ). Finally, figures are produced with the GREG program in the GILDAS package 0. When 



intensity ratios are presented, we use units of K km s"^ instead of Jy beam"' km s"' . Here we 
present the results for Model 2 and 4 since Model 2 is simpler to interpret than Model 1 , which has 
a "torus," and Model 4 has significantly different results from other models. 

Model 2: Because of the radial differences in abundances and excitation conditions, the ap- 
parent peaks in the intensity of molecular emission can change depending on the species and the 
line considered. Transitions of CO 7=1-0, 7 = 2-1, and 7 = 3-2 peak around lOpc from the 
AGN core with velocity-integrated intensities of about ~ 3000- 3500 K km s"' (see Fig [13] up- 
per panels). In general, since higher transition lines have higher critical densities, their intensities 
decrease with increasing radius faster than lower transition lines. For example, we find that the 
CO(3-2)/(l -0) ratio is slightly higher than 1 — around 1.2 — within 30 pc of the AGN, but 
falls below unity farther away. The HCN 7=1-0 and 7 = 3-2 emission peaks have intensities 
of 2500 and 2200 K km s"' as it can be compared in the upper panels of Fig [15] The ratio of 
HCN(7 = 3 - 2)/(7 = 1 - 0) is highest at 100 pc away from the AGN with a value of 1 .4. This ratio 
is smaller (~0.8) closer to the AGN. Since the CN abundance is lower near the AGN, the peak of 
CN (25/2-1 3/2) is about lOOpc away from the AGN with an intensity of 90 K km s"' (see Fig[T5]left 
panel). The HC0'^(1 -0) emission is about a half of HCN(1 -0) on scales of a few to 10 pc, but 
the ratio HC0VHCN(2- 1) becomes close to unity with increasing radius as shown in Fig [16] The 
intensity of HC3N(10-9) is on the order of thousands of K km s~\ equivalent to the HCN emission 
intensity. Intensity ratios of HC3N(10-9) over HCN(l-O) can be only sl ightly lower than the order 



of unity in (U)LIRGs (e.g, ~ 0.78 in Arp 220 and 0.40 in NGC 4418 iLindberg et al.ll2011l) . but 
not in typical AGNs. The abundance may instead be suppressed because of the grain depletion. 
Alternatively, the bottom right panel of Fig [12] shows that even a moderate amount of cosmic -ray 
ionization can dissociate HC3N and also suppress its abundance. 

Model 4: Since the overall column density of CO is much lower than that of Model 2, the 
peak velocity-integrated intensity is also much lower. The peak intensity of CO(l-O), C0(2- 
1), and CO(3-2) are 200, 300, and 350 K km s"'. The intensity ratios of CO(3-2)/(l -0) and 
C0(2- 1)/(1 -0) are 1.5 and 2 respectively around the AGN. The velocity-integrated intensity of 
HCN(1 -0) is ~ 150 K km s"' within 30 pc, but it falls rapidly with increasing distance from the 
AGN, and is around 50 K km s"' near 50 pc. The intensity ratio of HCN(3 -2)/(l -0) is about 2 
around 10 pc from the AGN, but it decreases with increasing radius, and it is below unity when 
r ~ 50pc. Much lower intensities than those of CO and HCN are predicted for CN. The intensity 
of CN(l3/2-0i/2) is 20K km s"' at lOpc from the AGN, and drops to 5 K km s"' near 20 pc. 



5http://www.iram.fi-/IRAMFR/GILDAS 
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CN(25/2- I3/2) has about 1.5 times more intensity than CN(l3/2-0i/2) in around lOpc, although 
it drops faster further out from the AGN. The higher level transition CN(37/2-25/2) has a higher 
intensity of 80 K km s~\ but decreases to less than 1 K km s"' farther out. Similar intensities to 
HCN are found for HCO'^(3-2) and (1-0) with the same 7-number transition near lOpc from the 
AGN core, but their intensity drops to about a tenth of the HCN lines around 30 pc away from the 
AGN. Although the abundance of HC3N is on the order of 10"^, which is much smaller than that 
in Model 2, there is a small region of several pc in Model 4 where HC3N emission is seen with 
HC3N(10-9)/HCN(l-0)~0.1. 

In general, almost the entire circumnuclear disk (CND) of Model 4 is XDR-dominated, and 
has higher intensities for the higher 7-transition lines compared with Model 2. Just as the abun- 
dance ratio of CN/HCN is higher in Model 4 than other models, the intensity ratio of CN/HCN is 
also higher. 



4. Discussion 

4.1. Comparison with Observations: A Case of NGC 1068 - a Prototypical AGN 

Although our models have not been tuned to represent any individual galaxy or AGN disk, 
a general comparison of models with observations helps to understand the physics and chemistry, 
and can serve as a guide both to future modeling and observations. Here, we compare our models 
with observations of a prototypical AGN containing galaxy NGC 1068. As a reference, observed 
intensities in NGC 1068 are listed in Table [T2l It must be noted that some of the observations in 
this table have different beamsizes and thus cannot be compared directly. 



4.1.1. CO 



The CO observations by iHelfer & Blitzl (|1995|) show sp iral arms on kpc scales, and a nuclear 
"ring" at radius ~ 200 pc. A higher resolution image by ISchinnerer et al.l (|2000|) reveals two 
peaks, each ~ 100 pc away from the central black hole, the so-called "Eastern" and "Western" 
knots, with Eastern knot having higher peak intensities by a factor of 2 than the Eastern knot for 
J = 1-0,2-1, and 3-2 lines of CO as summarized in Table [121 If the disk is a regular Keplerian 
disk and the CO abundances do not vary strongly, then the high velocity components must come 
from regions close to the center of the CND, and the two emission peaks should be symmetric. 
However, the CO emis sion appears asymmetric, and off from the expectations of a simple disk. 
Schinnerer et al.l (|2000|) suggest a warp in the disk as a possible cause for the observed double 
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knots with a complicated structure. iGalliano et al.l (|2003|) constructed a disk model without a 
warp, which is irradiated by X-ray s with an angle offset from the perpendicular axis of the disk. 
According to iMaloney et al.l (|1996|) . the CO intensity peaks where the X-ray intensity per density 



is moderately high, with Hx/n ~ 10"^^- 10"^^ erg cm^ s"' (See Eq. IA2l for the definition of Hx)- 
If the peaks are indeed created at regions illuminated by the X-ray ionizati on cones, one of them 
should appear more obscured th an the other, as illustrated in Figure 9 of IGalliano et al.l (|2003|) . 
A hydrodynamic simulation by IWada & Tomisakal (l2005h without X-ray irradiation and with a 
constant value for the so-called X-factor, the conversion factor from CO intensity to H2 column 
density, also explains this asymmetric double peak with sufficiently large CO beams of ~ 1", 
although smaller beams at a fraction of an arcsecond are predicted to reveal ring-shaped emission 
in their model. 

In our models, we find that the CO fractional abundance is mostly constant with (5-7) x 10"^ 
throughout the entire molecular region, and higher molecular emission also seem to come from the 
regions with larger column densities. Our models 2 and 4 both show centrally-peaked emission, 
which does not match the observations of NGC 1068. Although some of the emission near the 
AGN core might be absorbed by the surrounding colder gas, it is likely that there is simply a lack 
of molecular gas in this r egion , perhaps caused by feedback from the AGN. The total gas mass 
presented in Hicks et al.l (|2009|) from their warm H2 observation within 60 pc from the AGN is 
10^ M0, which is about an order of magnitude lower than the case in Model 2. 

The ratio between high-7 and low-7 transitions — e.g., CO(3-2)/CO(l -0) — is higher in 
the observations than in Models 2 and 4, which could also be explained if the density is lower than 
our models since there will be smaller amount of cold gas with the lower density. In Model 2, 
this ratio does not exceed unity in units of K km s"^ and in Model 4 these ratios are only ~ 1.5 
around their peak location. As shown in Table [I2l the observed ratio of C0(3 - 2)/( 1 - 0) in NGC 
1068 is 7. Although Model 4 is closer to the observations in terms of the ratio of high-7 to low-7 
transitions, the intensities fall more rapidly than the observations with increasing radius. Model 
2, on the other hand, has too high an intensity near the AGN core, but matches the observations 
fairly well on larger scales. Thus the actual CND of NGC 1068 must have lower density than our 
models on pc scales, but similar structure to Model 2 on 10 pc scales. Of course, because we have 
not attempted to fit the observations of NGC 1068, only a qualitative comparison should be made. 
We discuss other molecules below, highlighting qualitative and interesting differences between the 
models and the observations as a guide to future work. 
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4.1.2. HCN&HCN/CO 



Similar to CO, the observed HCN emission also has Western and Eastern knots, with stronger 
emission fr om the former . Alth ough detailed spatial information for HCN within the CND is 
not known, iTacconi et al.1 (119941) observed the HCN(1 -0) at a similar level of intensity to CO 
(1 -0). For 7 = 3-2, the HCN/CO int ensity ratio in uni ts of K km s"' is around ~ 0.2, and 
increases to ~ 0.3 around the AGN core. lUsero et al.l (|2004l) estimated the column density ratio of 
N( liCN)/N(CO) to be (1.6- 2.0) x 10"^, which is close to the value from a more recent observation 



bv lKamenetzkv et alJ (1201 lb . who estimated that A^(HCN)/A^(CO) = 1 .2 x 10"^ iKrips et all (1201 lb 
estimated slightly lower column density ratio of A'^(HCN)/A'^(C0)=3x lO"'*. 

These estimates assume a single component of homogeneous density and temperature, and 
different physical conditions within the CND cannot be taken into consideration. 

Our results in Figure |9] show that the calculated A^(HCN)/A^(CO) column density ratio can be 
more than one order of magnitude higher than the observed value for r < 20-30 pc, but it is 
within an order of magnitude to the observed ratio for larger r att = Across- The steady-state column 
density ratio of N(liCN)/N{CO) is at least an order of magnitude lower than the observed value for 
r > lOpc in Model 3 and 4. In Models 1 and 2, the ratio stays around the observed value because 
of the flared structures, as mentioned in Section 13.1.31 

In our models 2 and 4, and for transitions 7 = 3-2 and 7 = 1-0, there are regions where the 
intensity ratio of HCN/CO is larger than unity. Those are the regions where the X-ray irradiation 
is stronger than other regions. For 7 = 3-2, the region where the intensity ratio HCN/CO >1 is 
physically smaller, and this will lead to lower HCN/C0(3 -2) with respect to HCN/C0(1 -0) if 
observed with a large beamsize. The peaks of both molecules reside within 10 pc in our models, 
which is again different from NGC 1068. 



4.1.3. HCO+/HCN 



Similar to HCN and CO, both the Western and Eastern knots are seen for HC0'^(7 = 3-2 
and 4-3) (see Fig. 2 iKrips et al.ll201 11) . Unlike HCN or C O, HCO^ ernission does not seem to 
have a high intensity red-shifted component that traces a jet (IKrips et al.ll201l[) . The emission is 
likely to be originating in the XDR regions since it peaks sligh tly closer to the AGN core than the 



suggested that t he col umn density 



Kamenetzkv et a 



(1201 If) shows that 



CO and HCN peaks by a few tens of pc. lUsero et al.l (|20C 
ratio of MHCO+)/^(HCN) ~ 0.6-1.3 while analysis by 
A^(HCO"^)/A^(HCN) = 0.06, which is an order of magnitude lower. IKrips et all (|2011l) also obtain a 
value of column density ratio HCO+/HCN that is lower than unity, i.e., HCO+/HCN ~ 2 x 10"^ - 
0.1. In most of our models, HCOVHCN abundance ratios at t = ?rrn« are somewhere between 
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the two observed abundance ratios, but the steady-state abundance ratios at r > lOpc are greater 
than indicated by these observations (see Fig. [TT). Besides the non-steady-state abundances at 
t = ?cross5 higher metal abundances can also suppress HCO"^, causing a lower HCO'^/HCN column 
density ratio. If some fraction of the HCN emission is coming from the jet, as discussed above, 
the estimated HCO'^/HCN ratio in the XDR layer could in principle be higher than so far observed, 
and a high abundance ratio of HCO'^/HCN in our results in some models cannot be excluded. 

In the nuclear region (r < 90 pc) of the disk in Model 2, HCO+/HCN (1 - 0) < 1, and its value 
becomes close to unity for a larger radius. The value in Kohnd (l2005b . which is the average of 
the entire CND, is 0.5 for NGC 1068. When our results are seen with a larger beamsize, the ratio 
becomes close to unity. 



4.1.4. CN, CN/CO, CN/HCN 



Based on the maps of iGarcfa-Burillo et all (120081) and iGarcfa-Burillo et all (l2010l) . CN also 
has two peaks, Westem and Eastern knots. The location of the CN/CO intensity peak is seen 
at about 0.5 arcsec (~ 36 pc) away from the AGN core while the CN/HCN velocity-integrated 
intensi ty peak is about 0. 1 at the AGN core and increases to 0.3 at a peak seen at about 1 arcsec (~ 
72 pc). lusero etal](l20o4 estimated the Af(CN)/A^(CO) column density ratio to be (1 .8- 10) x 10"^ 
from their CN(2-1) observation, but more recent observations estimate a lo wer overall column 
density ratio of A^(CN)/A^(CO) = (0.3 - 1 .3) x lO'^ dGarcia-Burillo et al.ll201ol )FI In our models, the 
column density ratio A^(CN)/A^(HCN) is in general higher at larger r in Models 1 and 2 at both 
times, while for A^(CN)/A^(CO) the dependence is far less pronounced, as can be seen in Figures |7] 
through [8l The higher intensity of CN with respect to HCN in the colder and outer part of the 
disk may come from the higher A^(CN)/A^(HCN) column density ratio due to the combination of 
non-steady-state chemistry (t ~ tdy„) and the effect of X-rays. Our Model 2 calculation also has a 
higher CN/HCN intensity ratio with increasing distance from the AGN while Model 4 has a higher 
ratio near the AGN core. When the region is an XDR, it must have a higher CN/HCN intensity 
ratio. A lower intensity ratio of CN/HCN may indicate that the X-rays may not be the only source 
of heating, but there might also be thermal heating (i.e., from coupling of gas with warm dust, or 
heating from shocks). 



^The CN fract ional abundance is listed in terms of CN/H2, so we use CO/H2= 8x10 ^, which was used in 
Usero et all (Eqq3), to obtain CN/CO. 
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4.1.5. CS 



The fractional abundance of CS that 



Usero et al.l (20041) estimated is 1 x 10 ^. This is an order 



of magnitude higher than the value in Martin et al.l (|2009l) . Both values are relatively high with 
respect to our models, except for Models 3 and 4 on 100 pc scales if we utilize the degree of 
sulfur d epletion onto the gr ains used standardly in the chemical models of dense clouds in our 
Galaxy (IGraedel et al.lll982h . In galactic PDRs, there is evid ence that sulfur depletio n is less than 



in regular cold dense clouds, making S/H~ (2-5) x 10 ^ (|Goicoechea et al 



2006 ). In diffuse 



clouds, there a r e claim s that the sulfur abundance is even higher than the solar value. For example, 
Lucas & Liszti (|2002|) show that the ratio of elemental abundance S/0 is four times higher in a 
diffuse cloud towards COph than that in a solar neighborhood. The observed abundance of CS 
cannot be achieved in our models unless there are shock waves to cause early-time chemistry (i.e., 
abundances at t = Across), there is less depletion of sulfur onto the grains than in our models, or the 
metallicity is even higher than the solar value. Super-solar metallicities are commonly inferred 
in other AGNs, an d it is likely that the nuclear environment in NGC 1068 is similarly metal-rich 
dpields et allboovl) . 



4.1.6. Carbon-Chain Molecules 



Carbon-chain m olecules HC ^N (iCostaglio la et al.J 120111) . C2H (ICostagliola et al.l 12011 



Nakajima et al.ll201 lb . and C-C3H2 (iNakajima et al 201 1.) have been detected in NGC 1068. Emis- 
sions from carbon-chain molecules are in general weaker than the aforementioned molecules, and 
there is little spatial information within the CND. Although HC3N was detected in NGC 1068, the 
line ratios of HC3N/HCN are around 0.2-0.4 in (U) LIRGs, which are much higher than the other 
types of galaxies, which have a ratio less than 0.05 (ICostagliola et al.ll201 1[) . The reason for the 
enhancement of HC3 N is unknown, but concentration of its emission in galactic nuclei is seen for 
the case of Mrk 231 dAalto et allboul) and NGC 4418 (Costagliola et al. A&A, submitted). In 
addition, the ALMA early science observation by Takano et al. shows HC3N in the CND, not in 
the starburst ring. Our results show the enhancement of HC3N at high temperature; nevertheless, it 
seems that HC3N is more susceptible to X-ray irradiation than C-C3H2, as can be seen by compar- 
ing results in Tables |4l |5l and |6l Although the central concentration of HC3N in the CND of NGC 
1068 is seen with the ALMA early science observation, the emission is much weaker than the our 
Model 2 predicts according to Takano et al. (in prep). 
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4.2. A Case for (U)LIRGs 



Luminous infrared galaxies (LIRGs) are star-forming galaxies with infrared luminosities of 
LiR > 10^^ Lq (Lir > IO^^Lq for ultra-luminous infrared galaxies (ULIRGs)). It is believed that high 
infrared luminosities are due to star formation, and t hat these galaxies are potentially an evol ution- 
ary phase in galaxie s that precedes the AGN phase ( Sanders et al. 1988 : Sanders & Mirabeli ri996: 



Hopkins et al.ll2006l) . In ge neral, they have a higher gas fraction, a lower metallicity, and a highe r 



star formation efficiency (jPownes & Solomonlll998[ iRupke et al.ll2008l : iGao & SolomonI 12004 ) 
In addition to an e nhancement of HC3N, elevated abundances of C2H2 are seen in (U)LIRGs by 
Lahuis et al.l (|2007r) . In our models with 2=1, some complex molecules such as C2H2 and HC3N 



are predicted to have elevated abundances in particular regions of the disk discussed in Section 
13.11 These molecules indicate the presence of a high-temperature medium with a relatively low 
ionization rate. At warm temperatures around 200- 300 K, the steady-state abundances of these 
molecules are low, but the peak abundances at earlier times are still enhanced. Observational stud- 
ies of C-C3H2 in (U)LIRGs might be useful and may give further information on the condition of 
the interstellar medium in these galaxies. A detailed discussion of the relation between molecular 
abundances and physical conditions in (U)LIRGs will be presented in a future paper. 



4.3. Variation in Physical Conditions 



There are many models for the accretion disks of galactic nuclei, and there are still many 
unknowns in some of the physical processes, such as the mass accretion rate onto the supermassive 
black hole, the cosmic -ray energy density, the star formation rate, and the disk thickness (i.e. the 
scale height h{r)). The scale height of the dust torus in NGC 1 068 is found to be surprisingly 
high ( its thickness is 2.1 pc over 3.4 pc radius (|Jaffe et al.l 120041) ) in the central few pc with no 
conclusively identified mechanism to support the disk vertically against its gravity. 



The model by lThompson et al.l (|2005r) . on which our calculation for midplane density is based, 
assumes a Toomre Q parameter of close to unity, while the disk hei ght is supported by radiation 



pressure on the dust from OB-stars. On the other hand, the m odel bylWada et al. (120091) is based on 
feedback from the supemovae and stellar winds described in Wada & NormanI ( 200 1 ) . This model 
does not assume 2 ~ 1, an d the midplane density w ithin r < 30pc is at least an order of magnitude 



lower than in the model of Thompson et al. ( 20051) . although it is difficult to compare a full hydro- 
dynamic model containing str ong time-dependen t density inhomogeneities to our simplified and 
smooth disks. In the model by IWada et al.l (2009), the r nean density at a few tens of pc is around 
10"^ cm"^ as opposed to the model by [Thompson et al.l (|2005|) where the mean midplane density 
varies from 10^-10^ cm"^. At the lower densities, the total gaseous column density through the 
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model is around 0.3, while our results have molecular fractions close to unity. In the IWada et al 



midplane is averaged to be a little less than 10^^ cm ^, and a significant amou nt of higher energy 
X-rays can still penetrate into the midplane. The mean fraction of H2 in the IWada et al.l (120091) 



(|2009l) model, the medium is i nhomogeneous, and the molecu les can exist only in higher-density 
"clumps". Another model, by ^ Schartmann et al. (2009, 2010 ). also contains an inhomogeneous 
medium with a lower midplane density of lO'^ cm"^ within 10 pc, and even lower farther out from 
the AGN core. The chemistry in this third model would also be strongly affected by the larger 
XDR. 

The physical conditions change depending on the evolutionary stage and the dominant pro- 
cesses at work inside the CND. Comparing our models with the observations of NGC 1068, the 
molecular density of our models on scales smaller than 10 pc appears much higher. The CND of 
1068 may thus simply have a lower gas fraction in this region, perhaps as a result of preceding 
AGN feedback and/or stellar feedback. If the density near the nucleus is low enough to produce 
a "hole" in the molecular disk, the location of the XDR will change, and the peak locations of 
each molecular line change. Our models with 2=1 (Models 1 or 2) fail to capture these time- 
dependent effects, and thus probably have much higher density than the CND of NGC 1068. The 
assumption that 2=1 may not hold in parts of the disk where it is, or has bee n, strongly dis- 
turbed by AGN activity or recent stellar feedback, as shown in ISani et al.l (l2012h. However, the 



assumption of Q ~ 1 appea rs to hold true in many of galactic disks (e.g., Martin & Kennicut^boOl 



Downes & Solomonll 19981) . Effects of a systematic and time-dependent variation in the Toomre Q 
paramter in more models must be studied in order to understand the chemical abundances and star 
formation condition in active galaxies. 



5. Summary and Future Work 

The chemistry in the molecular disk around active galactic nuclei can be affected by many 
physical components of the system, including X-rays from the AGN core, the disk density struc- 
ture, star formation activity and cosmic rays, and dynamical processes such as shocks, turbulence, 
rapid accretion, and jets. In this paper, we considered simple density structures for the molecular 
disk, with X-rays, cosmic-rays, and gas-grain collisions as the primary heating mechanisms for the 
gas. We calculated the time-dependent abundances in the gaseous disk, reporting the abundances 
at both the characteristic shock-crossing time at each disk annulus and in steady state, and we cou- 
pled this abundance determination to a self-consistent calculation of the gas temperature. Although 
shocks are considered in the sense that they dissociate molecules, we did not include shock heating 
and the interdependence between this heating and the chemistry in this paper, which involves more 
details. 
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Our models of AGN disks show that there is an XDR layer, where the molecules are irradiated 
by X-rays, and a midplane, where the gas is colder, and more fully shielded, and it is this reservoir 
of mass that may be accreted into the black hole. Depending on the disk density structure and 
X-ray illumination, the location of the peak abundances, emission, and line ratios — both between 
different rotational levels and different molecular species — change. For example, the abundance 
ratio of CN/CO can be enhanced in the XDR layer. In addition, the abundance of CN in the 
crossing time is enhanced with respect to the steady-state abundance, and the abundance ratio of 
HCO'^/CO can also be enhanced in XDR layer with little time dependence. Time dependence is 
also important for HCN; the abundance ratio of HCN/CO can be enhanced in a warm part of the 
midplane at t = Across- We also provided preliminary radiative transfer calculations in our model 
disks and abundance profiles in order to make a first comparison with observations. The more the 
fraction of XDR layer the model has, the higher the ratio of higher-7 line intensity over lower-7 
line intensities. The enhancement of abundance in CN and HCO"^ in the XDR layer also appear as 
higher intensities. Our goal in the future is to present the full modeled intensities, and to apply the 
calculations in this paper to many different disk models in order to fit the observations in specific 
galaxies. There are also parameters that need more careful consideration such as the dumpiness, 
and variation in metallic ity. Although existing observatio ns can provide high-resolution intensity 
maps of CN and HCN (|Garcia-Burillo et al.ll2008l |2010|) in 10-pc scale, even higher resolution 
of pc-/subpc-scale images are possible with ALMA. In addition, the ALMA early science also 
demonstrated its power to observe more complex molecules. Strong lines of HC3N(11-10) and 
(12-11) were observed. It will be interesting to observe higher-J lines of HC^N such as J = 25-24 
that was observed with high excitation temperature in NGC 4418 by ICostagliolal (|2012|) . Radial 
dependences of abundances for CN, HCN, HCO"^, sulfur-containing species, and carbon chain 
molecules can be used to constrain high/low density physical models, metallicity, and the effects 
of X-rays. In particular, with the radial abundance gradient, one could constrain the disk structure 
including density and temperature from the molecular observations. 
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A. Temperature Calculation 
A.l. Heating 



A. 1.1. X-ray heating 



When atoms and molecul es are ionized by X-r ays, fast electrons are produced, which can 
interact with ambient electrons. iMaloney et al.l(|1996|) have shown that this Coulo mb heating domi- 
nates heating when the electron fraction is high. The Coulomb rate is expressed as (iMeijerink & Spaans 
2OO5I) : 

T = r]nHx (Al) 



where n is the total hy drogen density, r? is th e heating efficiency as inlMeiierink & SpaansI (|2005|) 
based on the results by iDalgamo et al.l (|1999|) . and Hx is defined as (|Maloney et al.lll996() 



Hx= apa{E)F{E)dE 



(A2) 



where (7pa{E)F(E) is the photoabsorption cross section and F(E) is the local photon energy flux 
at energy E. The gas can also be heated when X-ray s ionize molecular hydrogen by a series of 
reactions. As described in lGlassgold & Langerl(|l973|) . the ionization heating rate is: 



7-1- H2 Hj + e 

H2 + H+^H+ + H 
m + e-^H2 + H. 



(A3) 

(A4) 
(A5) 



The last two reactions are exothermic, a nd a large fraction of the e nergy goes into heating. We 
follow the formula for the heating rate in IMeijerink & SpaansI (|2005|) : 



r= , . , . , ^ X 10 '\H,xn,nerg cm - s \ 



(A6) 



where he, kn, and kni are the total reaction rates of Equations IA51 IA4I and IA3I while Xe, xu, and 
xm are the ionization fraction, and the fractional abundances of H and H2, respectively. is the 
ionization rate per hydrogen molecule, and n is the density. 



A.l. 2. Cosmic-ray heating 



As noted earlier, cosmic -rays can ionize atoms and molecules, which can produce secondary 
electrons. These secondary electrons can excite molecular hydrogen, which can be used in heating. 
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With the amount of e nergy deposited per primary ionization bv lCravens & Dalgamol (Il978h and 



Glass gold & Langerl (|l973|) . the heating rate by cosmic rays is given by iTielens & HoUenbach 
Jl985l) to be 

r = 1 .5 X IQ-'^HiriH, erg cm"-' s'^ , (A7) 



where Cm is the cosmic -ray ionization rate per molecular hydrogen and nm is the H2 density. 



A.1.3. Dust-gas interaction 
When the dust temperature and the gas temperature are different, they c an be heated or cooled 



by eac h other through collisions. We utilize the heating/cooling rate given by lHoUenbach & McKee 
d 19891) : 



r= 1.2x lO-^^n^ (1000^) [~^J [l-0-8exp(-75^/r)](r-r,„,,)ergcm-3 s"', 

(A8) 

where n is again the gas density, T is the kinetic gas temperature, Tjust is the dust temperature, and 
amin is the minimum grain size. 



A.2. Cooling 



The vibrational and rotati onal cooling rates of H2, CO, and H2O are taken from 
Neufeld & KaufmanI (Il993f) and iNeufeld et al.l (Il995h . In addition, atomic line cooling by CII 
158 fim, 01 63/xm and 63 00 A is included. These cooling rates are calculated using the method 
described in lTielens ( 2005 ): 



n^A = 



gu/giexp(-Eui/kT) 



-AjnAuihu^il3(T), 



(A9) 



l + ncr/3(T)+gJgiexpi-Eui/kTy 

where n is the total number density, and ^1 are the statistical weights of the upper and lower 
levels, ^ui is the energy difference between the upper and lower levels, k is the Boltzman constant, 
Aj is the abundance of species j, Aui is the Einstein coefficient, and Uui is the frequency of the line 
of the transition from the upper state to the lower state. The critical density is defined as 

_ /3(r)A„, 



(AlO) 



where 0(t) is the escape probability at an optical depth r. We used the form of escape probability 
as follows (|Tielensll2005|) 

/?(r)= '-7<:"""' .(r<7) (All) 

4. DOT 
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4r 



In ^ 



1 0.5 ' 



(r > 7). 



(A12) 



The optical depth is defined by 



Snu^ b/Az [riugi 



nigu 



(A13) 



In Equation (IA13I) . b is the Doppler broadening, and Az is the distance from the surface. This 
quantity ofb/Az can be considered as a velocity gradient, and the value of Azf2 was used. 
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Table 1 . Initial Fractional Abundances With Respect to Total Hydrogen. 



Species Low Metal Solar 



Hz 


0.5 


0.5 


He 


0.14 


0.09 





1.76(-4) 


2.56(-4) 


N 


2.14(-5) 


7.6(-5) 


F 


2.0(-8) 


1.8(-8) 


CI 


3.0(-9) 


1.8(-7) 


C+ 


7.3(-5) 


1.2(-4) 


Fe+ 


3.0(-9) 


2.0(-7) 


Mg+ 


3.0(-9) 


2.4(-6) 


Na+ 


3.0(-9) 


2.0(-7) 


P+ 


3.0(-9) 


1.17(-7) 


S+ 


2.0(-8) 


1.5(-5) 


Si+ 


3.0(-9) 


1.7(-6) 



Table 2. List of Symbols. 



Symbols 


Meaning 


P 


gas mass density 


n 


gas number density 


n 


Keplerian rotation frequency 


Q 


Toomre Q parameter 


Mbh 


Black hole mass 


a 


velocity dispersion 


Vs 


sound speed 


K 


epicyclic frequency 


E 


surface density 


/g 


gas fraction 


h 


scale height 


i-AGN 


total AGN luminosity 


P* 


star formation rate per volume 


V 


star formation efficiency 


^dyn 


dynamical time 


Across 


crossing time (~ tdyn) 


£/cR 


Cosmic-ray energy density 




Cosmic-ray luminosity 


fpp 


proton-proton decay time 


^wind 


wind escape time scale 


CcR 


cosmic-ray ionization rate 


Cx 


X-ray ionization raet 


CTJ 


Thomson cross section 
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Table 3 . Radial Dependence of the Scale Height h / r. 



disk model number 


<lpc 


l-3pc 


3-lOOpc 


P(z)/p(0) 


1 


0.01(r/pc)2 


0.5 


0.01(r/pc) 


exp(-^), CQMO?" for l-3pc 


2 


0.01(r/pc)2 


0.01(r/pc) 


0.01(r/pc) 


2 

exp(-|^) 


3 


0.01(r/pc)2 


0.5 


0.5 


exp(-^) 


4 


0.01(r/pc)2 


0.5 


0.5 


CQM07forl-100pc 



"Partial power-law dependence of density on height by fitting IChang et al.l |2007|) (CMQ07) 



Table 4. Integrated Vertical Column Densities (cm ^) for Model \, Q=\ dXt = t^ 



Radius(pc) 1 3 6 16 40 100 

7;„id (K) 750 470 300 190 120 75 

fcross(yr) 3x10^ 8x10^ 2x10* IxlO' 2x10^ 5x10^ 



CO 


l.Oe+W 


6.1e-l-20 


2.6e-l-20 


3.1e-l-20 


2.4e-l-20 


1.3e-l-20 


c 


4.3e+16 


2.2e-l-19 


e.ie+n 


2.3e-l-18 


2.7e-l-19 


6.2e-i-19 


C02 


l.le+n 


6.0e-l-18 


2.0e-l-18 


1.5e-l-17 


7.9e-l-17 


2.9e-l-17 


OH 


3.9e+17 


1.2e+19 


3.2e+18 


6.4e+18 


4.2e+18 


4.1e+18 


H2O 


1.7e+18 


3.3e+20 


1.2e-l-19 


3.8e+18 


1.8e-l-18 


1.6e-l-18 


NH3 


1.6e+13 


9.2e+16 


4.8e+16 


7.1e-l-16 


4.0e-l-16 


3.3e-l-16 


HCN 


4.7e+16 


2.2e+18 


8.4e+17 


l.le+17 


2.4e+17 


5.9e+16 


HNC 


4.5e+15 


1.3e+17 


3.1e+17 


1.2e+17 


2.3e+17 


5.5e+16 


CN 


1.7e-l-14 


2.4e-l-16 


1.9e-l-16 


2.9e-l-16 


2.8e-l-16 


1.8e-l-16 


cs 


3.3e-l-10 


4.0e-l-13 


2.5e-l-13 


e.oe+u 


1.9e-l-14 


1.5e-l-15 


SO 


2.4e-l-15 


7.4e+16 


7.2e-l-16 


2.0e-l-16 


6.6e-l-15 


1.6e-l-15 


S02 


2.9e+14 


5.7e+16 


1.4e-l-16 


3.1e-l-15 


1.4e-i-15 


2.0e-l-14 


C2H2 


7.3e+16 


5.8e+17 


1.9e-l-18 


2.2e-l-17 


5.4e-i-14 


6.7e-i-14 


HC3N 


2.8e-l-17 


9.9e-l-18 


l.le-l-19 


l.le-l-18 


6.7e-i-15 


1.8e-i-15 


HCO+ 


2.2e-l-15 


1.9e-l-17 


l.Oe+n 


1.6e-l-17 


3.3e-i-16 


2.7e-l-16 


H3O+ 


i.ie+ie 


3.6e-l-17 


7.2e-l-16 


1.3e-l-17 


7.5e-i-16 


6.4e-i-16 


HCNH+ 


l.le+B 


1.4e-i-15 


1.5e-i-15 


9.8e-l-14 


2.5e-i-14 


2.1e+14 


NH+ 


8.0e-l-10 


2.6e-l-14 


l.le+U 


1.9e-l-14 


8.4e-i-13 


8.8e-i-13 


H2O+ 


1.3e-l-14 


2.3e-l-16 


1.2e-l-15 


3.3e-i-15 


7.7e-l-15 


6.1e-i-15 


OH+ 


8.1e-l-13 


1.2e-l-17 


3.7e-i-15 


1.2e-i-16 


l.Oe-l-16 


5.9e-i-15 


C2H 


5.8e+ll 


l.Oe+14 


8.0e+14 


2.4e+14 


3.7e+15 


1.6e+16 


C3H2 


l.Oe+n 


1.2e-l-17 


6.4e-l-17 


4.3e-l-17 


2.0e-l-14 


6.2e-i-14 


NO 


1.7e-l-16 


2.8e-l-18 


1.6e-l-18 


2.6e-l-18 


4.2e-l-17 


4.1e-i-17 




1.5e-l-16 


3.9e-l-17 


2.4e-l-17 


3.5e-l-17 


3.1e-l-17 


5.6e-l-17 
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Table 5. Integrated Vertical Column Densities (cm ^) for Model 2, Q =1 at t = f( 



Radius(pc) 


1 


3 


6 


16 


40 


100 


Tmii (K) 


750 


470 


300 


190 


120 


75 


tcross (yr) 




Owl a3 
oX 10 


Owl n4 
2X lU 


1 X lU 


2X 10^ 


c w 1 n5 
JXIU 


CO 


4.5e+20 


3.4e+20 


2.2e+20 


1.8e+20 


1.8e+20 


1.9e+20 


c 


6.9e+17 


2.7e+18 


1.5e+19 


2.4e+19 


1.5e+19 


1.9e+19 


C02 


5.6e+18 


4.4e+18 


5.3e+17 


1.2e+18 


4.7e+17 


2.5e+17 


OH 


4.9e+18 


4.4e+18 


l.Oe+18 


1.9e+18 


3.0e+19 


2.9e+19 


H2O 


7.0e+20 


5.1e+19 


1.6e+20 


l.le+20 


6.4e+19 


1.3e+20 


NH3 


l.le+16 


9.7e+15 


4.0e+16 


1.5e+17 


l.le+17 


7.4e+16 


HCN 


4.7e+18 


3.5e+18 


1.6e+18 


5.5e+17 


2.1e+17 


l.le+17 


HNC 


1.2e+16 


2.3e+17 


5.1e+17 


2.9e+17 


8.4e+16 


2.5e+16 


CN 


1.3e+15 


2.3e+15 


3.2e+14 


3.1e+15 


4.4e+15 


4.6e+15 


CS 


2.2e+14 


3.7e+ll 


2.5e+13 


1.7e+14 


1.3e+15 


5.5e+14 


so 


6.7e+16 


6.0e+16 


2.6e+16 


1.2e+16 


7.8e+15 


4.4e+15 


SO2 


4.0e+16 


4.0e+16 


3.1e+15 


1.7e+15 


1.8e+15 


6.9e+14 


C2H2 


2.3e+18 


2.0e+17 


1.3e+18 


9.1e+15 


4.6e+14 


6.7e+14 


HC3N 


9.3e+18 


1.6e+19 


9.7e+18 


6.2e+15 


2.0e+15 


4.7e+14 


HCO+ 


4.0e+16 


6.0e+15 


9.6e+16 


3.9e+16 


l.le+17 


4.8e+17 


H3O+ 


9.1e+17 


9.6e+16 


6.5e+17 


3.2e+17 


9.3e+17 


1.9e+18 


HCNH+ 


3.2e+16 


8.5e+13 


4.7e+15 


3.0e+15 


2.1e+15 


1.7e+15 


NHJ 


5.9e+14 


1.7e+12 


3.1e+14 


5.0e+14 


4.3e+14 


3.9e+14 


H2O+ 


4.6e+13 


8.9e+15 


1.5e+15 


1.5e+16 


1.7e+16 


1.4e+16 


0H+ 


9.9e+15 


7.6e+15 


7.0e+16 


5.0e+17 


9.1e+17 


6.3e+17 


C2H 


1.7e+13 


4.7e+13 


1.2e+14 


1.8e+14 


2.1e+14 


4.2e+14 


C3H2 


9.0e+15 


1.4e+17 


6.2e+17 


3.8e+14 


2.0e+14 


3.1e+14 


NO 


5.1e+18 


6.6e+16 


1.4e+18 


1.2e+18 


3.3e+18 


4.0e+18 


HJ 


2.4e+17 


3.0e+16 


1.7e+17 


6.3e+16 


3.7e+17 


1.8e+18 
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Table 6. Integrated Vertical Column Densities (cm ^) for Model 3, Q=5 att = 



Radius(pc) 


1 


3 


6 


16 


40 


100 


Tmii (K) 


750 


470 


300 


190 


120 


75 


tcross (yr) 




Owl a3 
oX 10 


2X lU 


1 X lU 


2X 10^ 


c w 1 n5 
JXIU 


CO 


7.5e+21 


2.0e+21 


5.8e+20 


1.5e+20 


2.5e+19 


7.9e+18 


c 


2.8e+20 


2.8e+19 


l.le+20 


1.2e+20 


8.2e+19 


3.5e+19 


C02 


3.1e+19 


1.7e+19 


4.6e+17 


9.3e+17 


1.5e+17 


2.5e+16 


OH 


3.5e+20 


3.5e+19 


5.0e+18 


7.6e+17 


4.9e+17 


5.6e+16 


H2O 


4.8e+21 


1.5e+20 


4.1e+18 


7.2e+17 


2.4e+17 


5.4e+16 


NH3 


4.0e+18 


4.5e+17 


2.5e+16 


2.3e+16 


1.2e+15 


7.7e+14 


HCN 


4.3e+19 


1.3e+19 


4.7e+17 


3.5e+17 


4.2e+16 


5.5e+15 


HNC 


6.5e+18 


3.2e+18 


3.8e+17 


3.5e+17 


4.1e+16 


5.1e+15 


CN 


1.5e+17 


9.3e+16 


3.3e+16 


1.5e+16 


8.0e+15 


8.2e+15 


CS 


l.le+15 


3.0e+14 


1.9e+15 


8.4e+15 


l.le+16 


3.6e+15 


so 


6.5e+17 


4.4e+17 


1.3e+15 


6.1e+14 


2.9e+12 


4.4e+12 


SO2 


4.5e+17 


3.8e+16 


1.5e+13 


1.9e+13 


6.0e+09 


1.3e+10 


C2H2 


2.1e+19 


2.1e+19 


2.3e+18 


6.6e+16 


l.le+15 


1.6e+15 


HC3N 


1.7e+19 


4.8e+19 


5.9e+18 


8.7e+15 


9.6e+14 


1.2e+14 


HCO+ 


6.1e+18 


9.0e+17 


2.0e+I6 


1.8e+16 


1.2e+15 


5.2e+14 


H3O+ 


2.0e+19 


l.le+18 


l.le+17 


1.5e+16 


9.3e+15 


1.3e+15 


HCNH+ 


5.3e+16 


8.0e+15 


1.2e+15 


3.4e+14 


9.5e+13 


3.3e+13 


NHJ 


1.4e+16 


l.le+15 


4.0e+13 


6.0e+13 


2.5e+12 


2.0e+12 


H2O+ 


1.9e+17 


7.3e+15 


1.4e+16 


6.3e+15 


1.6e+15 


4.7e+13 


0H+ 


7.6e+17 


9.0e+16 


7.1e+16 


2.8e+16 


1.6e+15 


4.8e+13 


C2H 


7.0e+14 


2.6e+15 


6.6e+16 


2.7e+16 


1.2e+16 


6.1e+15 


C3H2 


5.9e+18 


1.7e+I9 


1.5e+19 


9.6e+15 


5.9e+14 


8.6e+14 


NO 


9.6e+19 


l.le+19 


5.1e+16 


1.8e+17 


2.3e+15 


3.9e+15 


HJ 


9.9e+18 


3.1e+18 


2.0e+17 


1.2e+17 


8.2e+16 


2.1e+16 
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Table 7. Integrated Vertical Column Densities (cm ^) for Model 4, Q =5 at t = f( 



Radius(pc) 


1 


3 


6 


16 


40 


100 


Tmii (K) 


750 


470 


300 


190 


120 


75 


tcross (yr) 




Owl a3 
oX 10 


2X lU 


1 X lU 


2X 10^ 


c w 1 n5 
JXIU 


CO 


1.2e+16 


3.2e+19 


1.4e+19 


4.4e+18 


1.2e+18 


3.5e+17 


c 


1.5e+19 


1.3e+19 


9.2e+18 


8.5e+18 


4.8e+18 


2.2e+18 


C02 


2.6e+07 


9.3e+15 


3.8e+16 


2.8e+16 


6.4e+15 


1.2e+15 


OH 


3.4e+16 


2.7e+18 


6.5e+17 


3.le+17 


l.3e+17 


1.2e+16 


H2O 


7.8e+14 


l.le+18 


1.9e+17 


6.8e+16 


2.5e+16 


4.4e+15 


NH3 


6.8e+05 


6.2e+15 


4.3e+15 


1.2e+14 


3.2e+13 


1.5e+13 


HCN 


4.0e+09 


9.4e+15 


1.5e+16 


l.Oe+16 


1.8e+15 


2.2e+14 


HNC 


6.5e+07 


7.6e+15 


1.5e+16 


l.Oe+16 


1.8e+15 


2.1e+14 


CN 


2.9e+12 


l.le+16 


2.2e+15 


9.7e+14 


6.9e+14 


2.4e+14 


CS 


2.3e+05 


6.3e+12 


2.7e+14 


8.0e+14 


6.0e+14 


2.0e+14 


so 


1.6e+07 


3.5e+14 


1.8e+14 


6.2e+ll 


1.3e+ll 


8.7e+10 


SO2 


l.Oe+02 


3.2e+13 


1.5e+13 


l.le+09 


2.6e+08 


2.2e+08 


C2H2 


8.4e+05 


2.7e+13 


1.2e+15 


1.6e+14 


7.6e+13 


3.4e+13 


HC3N 


4.6e-05 


4.8e+10 


4.0e+14 


2.0e+14 


3.8e+13 


5.2e+12 


HCO+ 


3.3e+ll 


3.3e+16 


4.7e+15 


4.8e+14 


8.3e+13 


2.4e+13 


H3O+ 


6.6e+13 


5.4e+16 


8.9e+15 


3.7e+15 


1.5e+15 


1.8e+14 


HCNH+ 


l.le+07 


2.7e+14 


4.0e+13 


9.5e+12 


3.9e+12 


l.le+12 


NHJ 


4.3e+02 


2.4e+13 


1.3e+13 


2.6e+ll 


6.8e+10 


3.7e+10 


H2O+ 


2.4e+15 


4.1e+15 


4.9e+15 


3.1e+15 


2.6e+15 


3.1e+13 


0H+ 


8.9e+16 


7.2e+16 


8.9e+16 


5.1e+16 


2.2e+16 


3.3e+13 


C2H 


3.4e+08 


7.5e+13 


7.5e+14 


1.2e+15 


7.8e+14 


2.5e+14 


C3H2 


4.6e+00 


2.9e+12 


8.6e+13 


l.le+14 


3.0e+13 


1.7e+13 


NO 


1.4e+13 


5.2e+17 


7.9e+16 


1.2e+15 


2.2e+14 


5.9e+13 


HJ 


1.3e+13 


1.2e+17 


3.5e+16 


2.2e+16 


l.le+16 


4.0e+15 
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Table 8. Integrated Vertical Column Densities (cm ^) for Model 1, 2 =1 at steady-state. 



Radius(pc) 
Tmid (K) 


1 

750 




3 
470 




6 

300 


16 
190 




40 
120 


100 
75 


CO 


9.4e+18 


5, 


.2e+2D 


3, 


,7e+20 


3.3e+20 


2 


.8e+20 


2.0e+20 


c 


2.3e+16 


2, 


.4e+19 


5, 


.8e+17 


1.5e+18 


3 


.2e+18 


3.7e+19 


C02 


5.6e+16 


2, 


.3e+18 


2, 


.6e+18 


3.5e+17 


1 


.6e+17 


5.7e+16 


OH 


5.6e+17 


1, 


.6e+19 


2, 


.7e+18 


7.2e+18 


1 


.le+19 


6.6e+18 


H2O 


2.1e+19 


8. 


.8e+20 


2. 


.8e+20 


7.7e+18 


5, 


.3e+18 


2.4e+18 


NH3 


1.7e+17 


5, 


.9e+18 


7, 


.3e+16 


9.7e+16 


7 


.Oe+16 


3.5e+16 


HCN 


1.3e+17 


4, 


.2e+18 


3, 


.9e+16 


4.0e+16 


4 


.4e+16 


1.2e+16 


HNC 


4.1e+16 


1, 


.5e+18 


1, 


.5e+16 


2.9e+16 


3 


.le+16 


9.0e+15 


CN 


3.0e+13 


1, 


,3e+16 


1, 


,le+16 


2.4e+16 


4 


.4e+16 


3.4e+16 


cs 


3.9e+14 


1 


,2e+16 


1. 


.8e+13 


4.6e+13 


5 


.4e+13 


6.3e+13 


so 


6.9e+14 


4, 


.le+16 


3, 


.Oe+16 


2.6e+16 


1 


.7e+16 


8.4e+15 


S02 


4.3e+14 


2, 


.3e+16 


5, 


.7e+16 


2.7e+16 


9 


.4e+15 


2.7e+15 


C2H2 


7.7e+17 


2, 


.6e+19 


1 


.3e+15 


8.6e+14 


4 


.Oe+14 


9.5e+14 


HC3N 


2.8e+15 


1 


.le+17 


1 


.3e+12 


8.()e+ll 


2 


.4e+ll 


1.3e+ll 


HCO+ 


4.0e+15 


2, 


.2e+17 


1 


.3e+17 


2.0e+17 


1 


.7e+17 


5.6e+16 


H3O+ 


2.7e+16 


8. 


.Oe+17 


1, 


.2e+17 


1.9e+17 


2 


.2e+17 


1. le+17 


HCNH+ 


2.1e+13 


2, 


,le+15 


7, 


.8e+14 


l.Oe+15 


1 


.2e+15 


2.8e+14 


NHJ 


2.4e+12 


5, 


,2e+14 


1, 


.7e+14 


2.3e+14 


1 


.7e+14 


7.8e+13 


H2O+ 


1.4e+14 


2, 


,2e+16 


2, 


,8e+14 


3.2e+15 


1 


.9e+15 


7.1e+15 


0H+ 


4.7e+13 


1 


.le+17 


2, 


.7e+15 


6.7e+15 


1 


.6e+15 


7.5e+15 


C2H 


5.6e+10 


4, 


.Oe+13 


5, 


.le+13 


1.6e+14 


3 


.4e+15 


1.5e+16 


C3H2 


1.3e+16 


2, 


.7e+17 


4, 


.6e+13 


4.2e+13 


1 


.5e+13 


9.8e+13 


NO 


5.7e+16 


3. 


.5e+18 


1. 


.8e+18 


3.1e+18 


2 


.6e+18 


6.3e+17 




6.9e+15 


3, 


.5e+17 


2, 


.3e+17 


3.5e+17 


1 


.le+18 


9.6e+17 
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Table 9. Integrated Vertical Column Densities (cm ^) for Model 2, 2 =1 at steady-state. 



Radius(pc) 


1 


3 


6 


16 


40 


100 


Tmid (K) 


750 


470 


300 


190 


120 


75 


CO 


5.5e+20 


3.4e+20 


3.6e+20 


3.0e+20 


2.7e+20 


2.0e+20 


c 


3.0e+18 


3.8e+17 


1.5e+18 


6.1e+18 


9.2e+18 


3.0e+19 


C02 


1.5e+17 


5.6e+18 


2.()e+18 


2.4e+17 


l.le+17 


3.1e+16 


OH 


1.4e+18 


1.2e+19 


6.8e+18 


1.2e+19 


1.2e+19 


l.le+19 


H2O 


7.9e+20 


7.1e+20 


2.6e+20 


1.2e+19 


5.7e+18 


3.6e+18 


NH3 


1.5e+19 


3.5e+18 


1.2e+17 


7.4e+16 


4.5e+16 


1.8e+16 


HCN 


4.1e+18 


5.2e+18 


l.le+17 


3.9e+16 


3.4e+16 


9.7e+15 


HNC 


1.2e+17 


1.5e+18 


1.7e+16 


1.7e+16 


2.2e+16 


6.6e+15 


CN 


6.8e+14 


8.2e+14 


9.4e+15 


2.2e+16 


3.8e+16 


3.0e+16 


cs 


1.3e+15 


5.6e+15 


3.1e+13 


4.1e+13 


3.4e+13 


4.8e+13 


so 


6.4e+15 


2.5e+16 


1.5e+16 


1.7e+16 


l.Oe+16 


4.7e+15 


S02 


9.9e+14 


2.6e+16 


4.4e+16 


1.6e+16 


5.5e+15 


1.5e+15 


C2H2 


1.3e+18 


1.6e+19 


8.2e+14 


4.9e+14 


3.0e+14 


9.4e+14 


HC3N 


l.le+17 


4.4e+16 


9.3e+ll 


3.9e+ll 


l.le+11 


l.Oe+U 


HCO+ 


l.()e+17 


1.5e+17 


2.2e+17 


2.()e+17 


1.4e+17 


5.1e+16 


H3O+ 


7.8e+17 


9.9e+17 


3.2e+17 


3.3e+17 


2.6e+17 


1.9e+17 


HCNH+ 


6.9e+15 


1.5e+15 


1.6e+15 


9.6e+14 


9.2e+14 


2.4e+14 


NHJ 


7.7e+14 


2.7e+14 


4.3e+14 


1.9e+14 


1.2e+14 


3.8e+13 


H2O+ 


3.8e+14 


1.7e+15 


1.2e+16 


1.3e+16 


1.5e+16 


1.8e+16 


0H+ 


2.0e+16 


2.8e+14 


l.le+16 


5.4e+16 


l.le+17 


l.le+17 


C2H 


1.8e+12 


9.0e+ll 


2.5e+13 


9.7e+13 


2.2e+14 


1.2e+16 


C3H2 


8.3e+15 


2.1e+17 


3.8e+13 


1.9e+13 


1.2e+13 


1.2e+14 


NO 


1.9e+18 


2.4e+18 


3.2e+18 


2.8e+18 


1.9e+18 


3.8e+17 




1.7e+17 


2.9e+17 


2.8e+17 


5.4e+17 


5.9e+17 


l.le+18 
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Table 10. Integrated Vertical Column Densities (cm ^) for Model 3, Q -5 at steady-state. 



Radius(pc) 


1 


3 


6 


16 


40 


100 


Tmid (K) 


750 


470 


300 


190 


120 


75 


CO 


7.8e+21 


2.7e+21 


9.3e+20 


3.4e+20 


1.2e+20 


4.8e+19 


c 


5.1e+19 


2.2e+19 


2.9e+18 


1.8e+18 


3.4e+18 


4.9e+17 


C02 


4.()e+18 


9.()e+18 


1.6e+18 


2.9e+17 


7.2e+16 


3.2e+16 


OH 


l.le+21 


4.4e+19 


8.3e+18 


2.3e+18 


1.5e+17 


2.1e+17 


H2O 


l.Oe+22 


2.9e+21 


1.6e+19 


1.7e+18 


3.2e+17 


2.3e+17 


NH3 


6.9e+19 


1.3e+18 


2.2e+17 


9.5e+16 


3.7e+16 


2.5e+16 


HCN 


3.9e+20 


l.Oe+19 


5.2e+16 


6.3e+15 


1.9e+15 


2.3e+15 


HNC 


6.5e+18 


5.4e+18 


3.9e+16 


5.4e+15 


1.6e+15 


2.1e+15 


CN 


2.8e+17 


2.8e+16 


2.9e+16 


4.6e+15 


3.5e+15 


6.0e+14 


cs 


7.7e+16 


1.9e+16 


6.9e+13 


3.8e+13 


2.0e+13 


3.4e+13 


so 


2.0e+17 


2.4e+17 


8.4e+16 


4.4e+16 


2.1e+16 


8.3e+15 


S02 


3.1e+16 


2.1e+17 


1.2e+17 


3.1e+16 


7.9e+15 


1.9e+15 


C2H2 


1.2e+20 


1.5e+19 


1.6e+15 


6.2e+14 


5.0e+14 


4.8e+14 


HC3N 


6.5e+17 


l.le+17 


l.le+12 


9.()e+l() 


6.8e+l() 


l.le+11 


HCO+ 


8.9e+18 


9.3e+17 


2.8e+17 


4.9e+16 


7.6e+15 


7.2e+15 


H3O+ 


4.4e+19 


2.7e+18 


2.1e+17 


5.0e+16 


3.7e+15 


4.6e+15 


HCNH+ 


l.le+18 


4.9e+15 


1.3e+15 


1.2e+14 


2.4e+13 


3.3e+13 


NHJ 


6.4e+16 


2.0e+15 


3.8e+14 


1.2e+14 


3.4e+13 


4.4e+13 


H2O+ 


2.0e+17 


9.7e+15 


1.2e+16 


1.7e+15 


7.0e+14 


2.0e+14 


0H+ 


2.8e+17 


9.9e+16 


1.7e+16 


2.1e+16 


3.3e+15 


1.9e+14 


C2H 


2.6e+18 


2.1e+13 


1.9e+14 


2.8e+13 


2.5e+14 


1.4e+13 


C3H2 


1.2e+19 


3.8e+17 


5.8e+13 


1.8e+12 


1.3e+12 


9.8e+ll 


NO 


8.0e+19 


1.3e+19 


4.4e+18 


7.9e+17 


8.6e+16 


1.4e+17 




2.7e+20 


2.0e+18 


4.1e+17 


8.2e+16 


2.5e+16 


1.2e+16 
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Table 11. Integrated Vertical Column Densities (cm ^) for Model 4, Q =5 at steady-state. 



Raclius(pc) 


1 


3 


6 


16 


40 


100 


Tmiii (K) 


750 


470 


300 


190 


120 


75 


CO 


1.2e+16 


3.3e+19 


2.6e+19 


1.5e+19 


6.4e+18 


2,7e+18 


c 


1.5e+19 


1.4e+19 


3.0e+18 


1.2e+18 


l.Oe+18 


1.4e+17 


C02 


2.6e+07 


1.2e+16 


2.4e+16 


l.le+16 


3.2e+15 


1.9e+15 


OH 


3.4e+16 


1.8e+18 


6.3e+17 


1.3e+17 


2.9e+16 


6.9e+15 


H2O 


7.8e+14 


3.1e+18 


3.2e+17 


5.6e+16 


1 .8e+16 


l.le+16 


NH3 


6.8e+{)5 


7.8e+15 


6.5e+15 


3.8e+15 


1.9e+15 


1.3e+15 


HCN 


4.()e+()9 


1.5e+16 


1.2e+15 


3.5e+14 


6.4e+13 


l.Oe+14 


HNC 


6.5e+07 


2.9e+15 


8.3e+14 


2.6e+14 


6.0e+13 


9.4e+13 


CN 


2.9e+12 


4.9e+15 


1.5e+15 


6.6e+14 


8.3e+13 


1.3e+13 


cs 


2.3e+05 


3.6e+12 


3.3e+12 


1.3e+12 


1.4e+12 


2.3e+12 


SO 


1.6e+07 


3.4e+14 


2.4e+15 


2.0e+15 


1.2e+15 


4.8e+14 


S02 


l.Oe+02 


4.0e+13 


1.7e+15 


1.3e+15 


4.2e+14 


l.Oe+14 


C2H2 


8.4e+05 


1.9e+14 


3.0e+13 


3.2e+13 


3.6e+13 


3.6e+13 


HC3N 


4.6e-05 


1.2e+ll 


l.le+10 


4.1e+09 


4.6e+09 


8.3e+09 


HCO+ 


3.3e+ll 


4.1e+16 


8.2e+15 


1.5e+15 


2.6e+14 


2.8e+14 


H3O+ 


6.6e+13 


5.9e+16 


1.2e+16 


1.5e+15 


2.9e+14 


1.8e+I4 


HCNH+ 


l.le+07 


2.8e+14 


2.8e+13 


7.0e+12 


4.3e+ll 


l.le+12 


NH+ 


4.3e+02 


3.5e+13 


1.3e+13 


3.8e+12 


1.4e+12 


2.0e+12 


H2O+ 


2.4e+15 


3.6e+15 


2.4e+15 


l.le+15 


3.0e+14 


6.0e+13 


0H+ 


8.9e+16 


6.8e+16 


5.0e+16 


3.9e+16 


l.le+16 


1.2e+14 


C2H 


3.4e+08 


1.5e+13 


6.1e+12 


4.9e+12 


8.7e+13 


6.5e+12 


C3H2 


4.6e+00 


1.2e+13 


5.0e+ll 


6.6e+10 


6.6e+ll 


8.3e+10 


NO 


1.4e+13 


4.8e+17 


1.4e+17 


1.9e+16 


2.7e+15 


4.3e+15 




1.3e+13 


8.4e+16 


1.5e+16 


6.5e+15 


4.0e+14 


4.2e+14 



Table 12. Observed Intensities in NGC 1068. Eastern knot is denoted as (E), and the western 

knot is denoted as (W). 



Species 


Lines 


Vel. Int. Peak Intensity 


Vel. Int. Peak Intensity 


Spatially Int. I 


Beam size 


Reference 






(Jy/beam km/s) 


(Kkm/s) 


(Jy km/s) 






i^CO 


1-0 


40(E), 19(W) 


760(E) 


90(E), 40 (W) 


1".0 X 0".8 


Krips etal. (2011) 




2-1 


70(E), 30(W) 


1980(E) 


290(E), 180(W) 


1".0 X 0".8 


Kripsetal. (2011) 




3-2 


470(E), 270(W) 


5900(E) 


1330(E), 720(W) 


1".0 X 0".8 


Krips etal. (2011) 


HCN 


1- 

2- 1 




121(E), 91(W) 


19(E),15(W) 


5" 


Tacconi et al.(1994) 




3-2 


51(E), 32(W) 


1100(E) 


110(E), 70(W) 


0".53x0".46 


Krips etal. (2011) 




4-3 




13.9 


220 


14" 


P.-B. et al. (2009) 


HCO+ 


1-0 














2-1 














3-2 


28 (E),14(W) 


590(E) 


52(E),40(W) 


1".0 X 0".8 


Kripsetal. (2011) 




4-3 


27(E) 




98(E) 


1".0 X 0".8 


Kripsetal. (2011) 


CN 


1-0 














2-1 


10(E), 3.6(W) 


200(E), 80(W) 









3-2 
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0.0 0.5 1.0 1.5 2.0 0.0 0.5 1.0 1.5 2.0 

lo5„(,odms(pt)) lo5,.(,odms(p=)) 



Fig. 1. — The radial and height dependence of density (color scale) and logarithms of column 
densities from the X-ray source \ogioN(cm~^) (contours). The radius, height, and the densities are 
all shown in logarithmic scale. Upper left panel: Model 1, Q=l; Upper right panel: Model 2, Q 
= 1; Lower left panel: Model 3, Q = 5 ; Lower right panel: Model A,Q = 5. 
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lo5.(,odi„s(p=)) Io,„(,adi„s(p0> 



Fig. 2. — The radial and height dependence of temperature r(K) (color scale) and logarithms of 
column densities from the X-ray source logioA^(cm"^) (contours). The radius, height, and temper- 
atures are all shown in logarithmic scale. Upper left panel: Model 1, Q =1; Upper right panel: 
Model 2, 2 = 1; Lower left panel: Model 3, 2 = 5; Lower right panel: Model 4, 2 = 5. 



Fig. 3. — The radial and height dependence of the ratio of the X-ray ionization rate over total hy- 
drogen density, Cx/^ (s"'cm^). The radius, height, and C/« are all shown in logarithmic scale. The 
XDR referred to in the text is shown in orange, red, and purple regions. The high-temperature 
region is within ~ 10 pc. Upper left panel: Model 1, Q=\\ Upper right panel: Model 2, 2 = 1; 
Lower left panel: Model 3,2 = 5; Lower right panel: Model A, Q = 5. 
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Fig. 4. — The radial and height dependence of fractional abundances ait = ta-oss in the Model 1, Q 
= 1 disk. Upper left panel: HCN; Upper right panel: CN; Middle left panel: HCO"^; Middle right 
panel: HC3N; Lower left panel: CS; Lower right panel: OVt. Note that the color scale is different 
forCS. 



-49- 




Fig. 5. — The radial and height dependence of fractional abundances at steady- state in the Model 
1, 2 = 1 disk. Upper left panel: HCN. Upper right panel: CH. Middle left panel: HCO"^. Middle 
right panel: HC3N. Lower left panel: CS. Lower right panel: 0H+. Note that the color scale is 
different for CS. 
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Fig. 6. — The radial and height dependence of fractional abundances ait = ta-oss in the Model 4, Q 
= 5 disk. Upper left panel: HCN; Upper right panel: CN; Middle left panel: HCO"^; Middle right 
panel: HC3N; Lower left panel: CS; Lower right panel: OVt. Note that the color scale is different 
forCS. 
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N(CN)/N(HCN) 



N(CN)/N(HCN) 





Fig. 7. — The ratio of CN/HCN vertical column densities vs radius. Left panel: at t = tr 



Right p anel: at steady-state. GB/U refers to the r atio o f the CN abundance by iGarcfa-Burillo et al. 



(120 101) over the HCN abundance by lUsero et al.l (|2004|) . U04 refers to the abundance ratio inferred 
bv lUseroetaP (12004 . 




Fig. 8. — The ratio of CN/CO vertical column densities vs radius. Left panel: at t = trm, , - Righ t 
panel: at steady-state. GBIO refers to the observed abundanc e by iGarcfa-Burillo et al.l (|2010|) . 
U04 refers to the abundance ratio inferred by lUsero et al.l(|2004|) . 
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N(HCN)/N(CO) 



N(HCN)/N(CO) 
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0.01 
0.001 
0.0001 
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100 
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Fig. 9. — The ratio of HCN/CO vertical column densities vs radius. L eft panel: att = trm../, Ri^h t 



panel: at steady-state. Kl 1 refers to the obser ved abundance of H CN by lKamenetzky et al.l (|201lh 
U04 refers to the abundance ratio inferred by lUsero et al.l(|2004l) . 



N(HCO+)/N(CO) 



N(HCO+)/N(CO) 



0.01 



0.001 - 



0.0001 : 



1e-05 : 



1e-06 : 



1e-07 



10 
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100 



0.01 



0.001 s 



0.0001 



1e-05 
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100 



Fig. 10. — The ratio of HCO'^/CO vertical column densities vs radi us. Left pane l : att = tcmss'. Right 
panel: at steady-state. U04 refers to the abundance ratio inferred b y lUsero et al.l(|2004h . Kl 1 refers 
to the observed abundance of HCO"^ by lKamenetzky et al.l (120 lib . 
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N(HCO+)/N(HCN) 
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Fig. 11. — The ratio of HCOVHCN vertical column densities vs radius. Left panel: ay_= tcross\ 
Right panel: at steady-state. U04 refers to the abundance ratio inferred bvlUsero et al.l(|2004() . Kl 1 
refers to the observed abundance ratio of HCO'^/HCN by lKamenetzky et a 
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Fig. 12. — The variation of assorted abundances vs radius with different cosmic-ray ionization 
rates caused by differing star formation efficiencies for the Model 1,2=1 disk. Upper left 
panel: HCN; Upper right panel: CN; Lower left panel HCO"^; Lower right panel: HC3N. 
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Fig. 13. — Velocity-integrated line intensities of CO emission are shown. Upper left panel: Model2 
CO(l-O). Contour levels are from 1500 to 3000 K km s"' in steps of 250 K km s"^; Upper right 
panel: Model 2 CO (3-2). Contour levels are from 1500 to 3500 K km s"^ in steps of 250 K km 
s"^; Lower left panel Model 4 CO (1-0). Contour levels are from 100 to 175 K km s"^ with 25 K 
km s"^ increments.; Lower right panel: Model 4 CO (3-2). Contour levels are from 100 to 350 K 
km s"' in steps of 25 K km s"'. 
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Fig. 14. — Velocity-integrated line intensities of HCN emission are shown. Upper left panel: 
Model2 HCN(l-O). Contour levels are from 1000 to 2500 K km s"' with each step of 250 K km 
s"^ Upper right panel: Model 2 HCN (3-2). Contour levels are from 1000 to 2000 K km s"' with 
each step of 250 K km s"'; Lower left panel Model 4 HCN (1-0). Contour levels are from 50 to 
125 K km s~' with each step of 25 K km s"'; Lower right panel: Model 4 HCN (3-2). Contour 
levels are from 100 to 225 K km s"' with each step of 25 K km s"^ 
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Fig. 15. — Velocity-integrated line intensities of CN emission are shown. Upper left panel: Model2 
CN(25/2-l3/2). Contour levels are 50, 65, and 80 K km s"^; Upper right panel: Model 4 CN (25/2- 
I3/2). Contour levels are 6, 8, 10, and 12 K km s"^ 
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Fig. 16. — Velocity-integrated line intensities of HCO"^ emission are shown. Upper left panel: 
Model2 HCO+(1-0). Contour levels are from 800 to 1200 K km s"' with each step of 100 K km 
s"^; Upper right panel: Model 2 HCO'^(3-2). Contour levels are from 1000 to 1750 K km s"^ with 
each step of 250 K km s"'; Lower left panel Model 4 HCO'^(l-O). Contour steps are 10, 20, and 30 
K km s"^; Lower right panel: Model 4 HCO+(3-2) Contour steps are 10, 20, 30, and 40 K km s"^ 



